text:
"""
Contains various miscellaneous useful functions.
"""
# load general packages and functions
import csv
from collections import namedtuple
from typing import Union, Tuple
from warnings import filterwarnings
import numpy as np
import matplotlib
matplotlib.use('Agg')
from matplotlib import pyplot as plt
import torch
import rdkit
from rdkit import RDLogger
from rdkit.Chem import MolToSmiles
from torch.utils.tensorboard import SummaryWriter

# load GraphINVENT-specific functions
from parameters.constants import constants

# defines the tensorboard writer
tb_writer = SummaryWriter(log_dir=constants.tensorboard_dir, flush_secs=10)

def get_feature_vector_indices() -> list:
    """
    Gets the indices of the different segments of the feature vector. The
    indices are analogous to the lengths of the various segments.

    Returns:
    -------
        idc (list) : Contains the indices of the different one-hot encoded
          segments used in the feature vector representations of nodes in
          `MolecularGraph`s. These segments are, in order, atom type, formal
          charge, number of implicit Hs, and chirality.
    """
    idc = [constants.n_atom_types, constants.n_formal_charge]

    # indices corresponding to implicit H's and chirality are optional (below)
    if not constants.use_explicit_H and not constants.ignore_H:
        idc.append(constants.n_imp_H)

    if constants.use_chirality:
        idc.append(constants.n_chirality)

    return np.cumsum(idc).tolist()

def get_last_epoch() -> str:
    """
    Gets previous training epoch by reading it from the "convergence.log" file.

    Returns:
    -------
        epoch_key (str) : A string that indicates the final epoch written to
          "convergence.log".
    """

    if constants.job_type != "fine-tune":
        # define the path to the file containing desired information
        convergence_path = constants.job_dir + "convergence.log"

        try:
            # epoch_key, lr, avg_loss
            epoch_key, _, _ = read_row(path=convergence_path,
                                       row=-1,
                                       col=(0, 1, 2))
        except ValueError:
            epoch_key = "Epoch 1"

        generation_epoch = constants.generation_epoch
        if constants.job_type == "generate":
            epoch_key = f"Epoch GEN{generation_epoch}"
        elif constants.job_type == "test":
            epoch_key = f"Epoch EVAL{generation_epoch}"
    else:
        # define the path to the file containing desired information
        convergence_path = constants.job_dir + "convergence.log"

        try:
            # epoch_key, lr, avg_loss
            epoch_key_tmp, _, _ = read_row(path=convergence_path,
                                           row=-1,
                                           col=(0, 1, 2))
            epoch_key_tmp       = epoch_key_tmp.split(" ")
            epoch_key = " ".join(epoch_key_tmp[:2])
        except (ValueError, IndexError) as e:
            epoch_key = "Step init"

    return epoch_key

def normalize_evaluation_metrics(property_histograms : dict,
                                 epoch_key : str) -> Tuple[list, ...]:
    """
    Normalizes histograms in `props_dict`, converts them to `list`s (from
    `torch.Tensor`s) and rounds the elements. This is done for clarity when
    saving the histograms to CSV.

    Args:
    ----
        property_histograms (dict) : Contains histograms of evaluation metrics
          of interest.
        epoch_key (str) : Indicates the training epoch.

    Returns:
    -------
        norm_n_nodes_hist (torch.Tensor) : Normalized histogram of the number of
          nodes per molecule.
        norm_atom_type_hist (torch.Tensor) : Normalized histogram of the atom
          types present in the molecules.
        norm_charge_hist (torch.Tensor) : Normalized histogram of the formal
          charges present in the molecules.
        norm_numh_hist (torch.Tensor) : Normalized histogram of the number of
          implicit hydrogens present in the molecules.
        norm_n_edges_hist (torch.Tensor) : Normalized histogram of the number of
          edges per node in the molecules.
        norm_edge_feature_hist (torch.Tensor) : Normalized histogram of the
          edge features (types of bonds) present in the molecules.
        norm_chirality_hist (torch.Tensor) : Normalized histogram of the chiral
          centers present in the molecules.
    """
    # compute histograms for non-optional features
    norm_n_nodes_hist = [
        round(i, 2) for i in
        norm(property_histograms[(epoch_key, "n_nodes_hist")]).tolist()
    ]
    norm_atom_type_hist = [
        round(i, 2) for i in
        norm(property_histograms[(epoch_key, "atom_type_hist")]).tolist()
    ]
    norm_charge_hist = [
        round(i, 2) for i in
        norm(property_histograms[(epoch_key, "formal_charge_hist")]).tolist()
    ]
    norm_n_edges_hist = [
        round(i, 2) for i in
        norm(property_histograms[(epoch_key, "n_edges_hist")]).tolist()
    ]
    norm_edge_feature_hist = [
        round(i, 2) for i in
        norm(property_histograms[(epoch_key, "edge_feature_hist")]).tolist()
    ]

    # compute histograms for optional features
    if not constants.use_explicit_H and not constants.ignore_H:
        norm_numh_hist = [
            round(i, 2) for i in
            norm(property_histograms[(epoch_key, "numh_hist")]).tolist()
        ]
    else:
        norm_numh_hist = [0] * len(constants.imp_H)

    if constants.use_chirality:
        norm_chirality_hist = [
            round(i, 2) for i in
            norm(property_histograms[(epoch_key, "chirality_hist")]).tolist()
        ]
    else:
        norm_chirality_hist = [1, 0, 0]

    return (norm_n_nodes_hist, norm_atom_type_hist, norm_charge_hist,
            norm_numh_hist, norm_n_edges_hist, norm_edge_feature_hist,
            norm_chirality_hist)

def get_restart_epoch() -> Union[int, str]:
    """
    Gets the restart epoch e.g. epoch for the last saved model state
    (`model_restart.pth`). Will simply return zero if called outside of a
    restart job.

    Returns:
    -------
        epoch (int or str) :
    """
    if constants.restart or constants.job_type == "test":
        # define path to output file containing info on last saved state
        generation_path = constants.job_dir + "generation.log"
        epoch           = "NA"
        row             = -1
        while not isinstance(epoch, int):
            epoch_key = read_row(path=generation_path, row=row, col=0)
            try:
                epoch = int(epoch_key[6:])
            except ValueError:
                epoch = "NA"
            row      -= 1
    elif constants.job_type == "fine-tune":
        epoch = constants.generation_epoch

    else:
        epoch = 0

    return epoch

def load_ts_properties(csv_path : str) -> dict:
    """
    Loads training set properties from CSV, specified by the `csv_path`, and
    returns them as a dictionary.

    Args:
    ----
        csv_path (str) : Path to CSV file from which to read the training set
          properties.

    Returns:
    -------
        properties (dict) : Contains training set properties.
    """
    print("* Loading training set properties.", flush=True)

    # read dictionary from CSV
    with open(csv_path, "r") as csv_file:
        reader   = csv.reader(csv_file, delimiter=";")
        csv_dict = dict(reader)

    # create `properties` from `csv_dict`, fix any bad filetypes
    properties = {}
    for key, value in csv_dict.items():

        # first determine if key is a tuple
        key = eval(key)
        if len(key) > 1:
            tuple_key = (str(key[0]), str(key[1]))
        else:
            tuple_key = key

        # then convert the values to the correct data type
        try:
            properties[tuple_key] = eval(value)
        except (SyntaxError, NameError):
            properties[tuple_key] = value

        # convert any `list`s to `torch.Tensor`s (for consistency)
        if isinstance(properties[tuple_key], list):
            properties[tuple_key] = torch.Tensor(properties[tuple_key])

    return properties

def norm(list_of_nums : list) -> list:
    """
    Normalizes input list of numbers.

    Args:
    ----
        list_of_nums (list) : Contains `float`s or `int`s.

    Returns:
    -------
        norm_list_of_nums (list) : Normalized version of input `list_of_nums`.
    """
    try:
        norm_list_of_nums = list_of_nums / sum(list_of_nums)
    except:  # occurs if divide by zero
        norm_list_of_nums = list_of_nums
    return norm_list_of_nums

def one_of_k_encoding(x : Union[str, int], allowable_set : list) -> 'generator':
    """
    Returns the one-of-k encoding of a value `x` having a range of possible
    values in `allowable_set`.

    Args:
    ----
        x (str, int) : Value to be one-hot encoded.
        allowable_set (list) : Contains all possible values to one-hot encode.

    Returns:
    -------
        one_hot_generator (generator) : One-hot encoding. A generator of `int`s.
    """
    if x not in set(allowable_set):  # use set for speedup over list
        raise Exception(
            f"Input {x} not in allowable set {allowable_set}. Add {x} to "
            f"allowable set in either a) `features.py` or b) your submission "
            f"script (`submit.py`) and run again."
        )
    one_hot_generator = (int(x == s) for s in allowable_set)
    return one_hot_generator

def properties_to_csv(prop_dict : dict, csv_filename : str,
                      epoch_key : str, append : bool=True) -> None:
    """
    Writes a CSV summarizing how training is going by comparing the properties
    of the generated structures during evaluation to the training set. Also
    writes the properties to an active tensorboard.

    Args:
    ----
        prop_dict (dict)   : Contains molecular properties.
        csv_filename (str) : Full path/filename to CSV file.
        epoch_key (str)    : For example, "Training set" or "Epoch {n}".
        append (bool)      : Indicates whether to append to the output file (if
                             the file exists) or start a new one. Default `True`.
    """
    # get all the relevant properties from the dictionary
    frac_valid  = prop_dict[(epoch_key, "fraction_valid")]
    avg_n_nodes = prop_dict[(epoch_key, "avg_n_nodes")]
    avg_n_edges = prop_dict[(epoch_key, "avg_n_edges")]
    frac_unique = prop_dict[(epoch_key, "fraction_unique")]

    # use the following properties if they exist
    try:
        run_time      = prop_dict[(epoch_key, "run_time")]
        frac_valid_pt = round(
            float(prop_dict[(epoch_key, "fraction_valid_properly_terminated")]),
            5
        )
        frac_pt       = round(
            float(prop_dict[(epoch_key, "fraction_properly_terminated")]),
            5
        )
    except KeyError:
        run_time      = "NA"
        frac_valid_pt = "NA"
        frac_pt       = "NA"

    (norm_n_nodes_hist,
     norm_atom_type_hist,
     norm_formal_charge_hist,
     norm_numh_hist,
     norm_n_edges_hist,
     norm_edge_feature_hist,
     norm_chirality_hist) = normalize_evaluation_metrics(prop_dict, epoch_key)

    if not append:
        # file does not exist yet, create it
        with open(csv_filename, "w") as output_file:
            # write the file header
            output_file.write(
                "set, fraction_valid, fraction_valid_pt, fraction_pt, run_time, "
                "avg_n_nodes, avg_n_edges, fraction_unique, atom_type_hist, "
                "formal_charge_hist, numh_hist, chirality_hist, "
                "n_nodes_hist, n_edges_hist, edge_feature_hist\n"
            )

    # append the properties of interest to the CSV file
    with open(csv_filename, "a") as output_file:
        output_file.write(
            f"{epoch_key}, {frac_valid:.3f}, {frac_valid_pt}, {frac_pt}, {run_time}, "
            f"{avg_n_nodes:.3f}, {avg_n_edges:.3f}, {frac_unique:.3f}, "
            f"{norm_atom_type_hist}, {norm_formal_charge_hist}, "
            f"{norm_numh_hist}, {norm_chirality_hist}, {norm_n_nodes_hist}, "
            f"{norm_n_edges_hist}, {norm_edge_feature_hist}\n"
        )

    # write to tensorboard
    try:
        epoch = int(epoch_key.split()[1])
    except:
        pass
    else:
        # scalars
        tb_writer.add_scalar("Evaluation/fraction_valid", frac_valid, epoch)
        tb_writer.add_scalar("Evaluation/fraction_valid_and_properly_term", frac_valid_pt, epoch)
        tb_writer.add_scalar("Evaluation/fraction_properly_terminated", frac_pt, epoch)
        tb_writer.add_scalar("Evaluation/avg_n_nodes", avg_n_nodes, epoch)
        tb_writer.add_scalar("Evaluation/fraction_unique", frac_unique, epoch)

def read_column(path : str, column : int) -> np.ndarray:
    """
    Reads a column from a CSV file. Returns column values as a `numpy.ndarray`.
    Removes missing values ("NA") before returning.

    Args:
    ----
        path (str)   : Path to CSV file.
        column (int) : Column to read from CSV file.

    Returns:
    -------
        data (np.ndarray) : Read column.
    """
    with open(path, "r") as csv_file:
        data = np.genfromtxt(csv_file,
                             dtype=None,
                             delimiter=",",
                             skip_header=1,
                             usecols=column,
                             missing_values="NA")
    data = np.array(data)
    data = data[~np.isnan(data)]  # exclude `nan`
    return data

def read_last_molecule_idx(restart_file_path : str) -> Tuple[int, int]:
    """
    Reads the index of the last preprocessed molecule from a file called
    "index.restart" located in the same directory as the data. Also returns the
    dataset size thus far.

    Args:
    ----
        restart_file_path (str) : Path to the index restart file.

    Returns:
    -------
        Tuple[int, int] : The first integer is the index of the last preprocessed
                          molecule and the second integer is the dataset size.
    """
    with open(restart_file_path + "index.restart", "r") as txt_file:
        last_molecule_idx = np.genfromtxt(txt_file, delimiter=",")
    return int(last_molecule_idx[0]), int(last_molecule_idx[1])

def read_row(path : str, row : int, col : int) -> np.ndarray:
    """
    Reads a row from CSV file. Returns it as a `numpy.ndarray`. Removes "NA"
    missing values from the column before returning.

    Args:
    ----
        path (str) : Path to CSV file.
        row (int)  : Row to read.
        col (int)  : Column to read.

    Returns:
        np.ndarray : Desired row from the file.
    """
    with open(path, "r") as csv_file:
        data = np.genfromtxt(csv_file,
                             dtype=str,
                             delimiter=",",
                             skip_header=1,
                             usecols=col)
    data = np.array(data)
    return data[:][row]

def suppress_warnings() -> None:
    """
    Suppresses unimportant warnings for a cleaner readout.
    """
    RDLogger.logger().setLevel(RDLogger.CRITICAL)
    filterwarnings(action="ignore", category=UserWarning)
    filterwarnings(action="ignore", category=FutureWarning)
    # could instead suppress ALL warnings with:
    # `filterwarnings(action="ignore")`
    # but choosing not to do this

def turn_off_empty_axes(n_plots_y : int, n_plots_x : int, ax : plt.axes) -> plt.axes:
    """
    Turns off empty axes in a `n_plots_y` by `n_plots_x` grid of plots.

    Args:
    ----
        n_plots_y (int) : Number of plots along the y-axis.
        n_plots_x (int) : Number of plots along the x-axis.
        ax (plt.axes)   : Matplotlib object containing grid of plots.

    Returns:
    -------
        plt.axes : The updated Matplotlib object.
    """
    for vi in range(n_plots_y):
        for vj in range(n_plots_x):
            # if nothing plotted on ax, it will contain `inf`
            # in axes lims, so clean up (turn off)
            if "inf" in str(ax[vi, vj].dataLim):
                ax[vi, vj].axis("off")
    return ax

def write_last_molecule_idx(last_molecule_idx : int, dataset_size : int,
                            restart_file_path : str) -> None:
    """
    Writes the index of the last preprocessed molecule and the current dataset
    size to a file.

    Args:
    ----
        last_molecules_idx (int) : Index of last preprocessed molecule.
        dataset_size (int) : The dataset size.
        restart_file_path (str) : Path indicating where to save indices (should
          be same directory as the dataset).
    """
    with open(restart_file_path + "index.restart", "w") as txt_file:
        txt_file.write(str(last_molecule_idx) + ", " + str(dataset_size))

def write_job_parameters(params : namedtuple) -> None:
    """
    Writes job parameters/hyperparameters to CSV.

    Args:
    ----
        params (namedtuple) : Contains job parameters and hyperparameters.
    """
    dict_path = params.job_dir + "params.csv"

    with open(dict_path, "w") as csv_file:
        writer = csv.writer(csv_file, delimiter=";")
        for key, value in enumerate(params._fields):
            writer.writerow([value, params[key]])

def write_preprocessing_parameters(params : namedtuple) -> None:
    """
    Writes job parameters/hyperparameters in `params` (`namedtuple`) to
    CSV, so that parameters used during preprocessing can be referenced later.

    Args:
    ----
        params (namedtuple) : Contains job parameters and hyperparameters.
    """
    dict_path = params.dataset_dir + "preprocessing_params.csv"
    keys_to_write = ["atom_types",
                     "formal_charge",
                     "imp_H",
                     "chirality",
                     "group_size",
                     "max_n_nodes",
                     "use_aromatic_bonds",
                     "use_chirality",
                     "use_explicit_H",
                     "ignore_H"]

    with open(dict_path, "w") as csv_file:
        writer = csv.writer(csv_file, delimiter=";")
        for key, value in enumerate(params._fields):
            if value in keys_to_write:
                writer.writerow([value, params[key]])

def write_graphs_to_smi(smi_filename : str,
                        molecular_graphs_list : list,
                        write : bool=False) -> \
                        Tuple[float, torch.Tensor, torch.Tensor]:
    """
    Calculates the validity and uniqueness of input molecular graphs. Then,
    (optional) writes the input molecular graphs a SMILES file.

    Args:
    ----
        smi_filename (str)           : Full path/filename to output SMILES file.
        molecular_graphs_list (list) : Contains molecular graphs.
        write (bool, optional)       : If True, writes the input molecular graphs
                                       to SMILES. If False, skips writing the SMILES,
                                       but calculates the rest of the metrics
                                       (validity, uniqueness) anyways.

    Returns:
    -------
        fraction_valid (float)           : The fraction of molecular graphs which
                                           are valid molecules.
        validity_tensor (torch.Tensor)   : A binary vector indicating the chemical
                                           validity of the input graphs using either
                                           a 1 (valid) or 0 (invalid).
        uniqueness_tensor (torch.Tensor) : A binary vector indicating the uniqueness
                                           of the input graphs using either a 1
                                           (unique or first duplicate instance)
                                           or 0 (duplicate).
    """
    validity_tensor   = torch.zeros(len(molecular_graphs_list),
                                    device=constants.device)
    uniqueness_tensor = torch.ones(len(molecular_graphs_list),
                                   device=constants.device)
    smiles            = []

    with open(smi_filename, "w") as smi_file:

        if write:
            smi_writer = rdkit.Chem.rdmolfiles.SmilesWriter(smi_file)

        for idx, molecular_graph in enumerate(molecular_graphs_list):

            mol = molecular_graph.get_molecule()
            try:
                mol.UpdatePropertyCache(strict=False)
                rdkit.Chem.SanitizeMol(mol)
                if write:
                    smi_writer.write(mol)
                current_smiles = MolToSmiles(mol)
                validity_tensor[idx] = 1
                if current_smiles in smiles:
                    uniqueness_tensor[idx] = 0
                smiles.append(current_smiles)
            except (ValueError, RuntimeError, AttributeError):
                # molecule cannot be written to file (likely contains unphysical
                # aromatic bond(s) or an unphysical valence), so put placeholder
                if write:
                    # `validity_tensor` remains 0
                    smi_writer.write(rdkit.Chem.MolFromSmiles("[Xe]"))

        if write:
            smi_writer.close()

    fraction_valid = torch.sum(validity_tensor, dim=0) / len(validity_tensor)

    return fraction_valid, validity_tensor, uniqueness_tensor

def write_training_status(epoch : Union[int, None]=None,
                          lr : Union[float, None]=None,
                          training_loss : Union[float, None]=None,
                          validation_loss : Union[float, None]=None,
                          score : Union[float, None]=None,
                          append : bool=True) -> None:
    """
    Writes the current epoch, loss, learning rate, and model score to CSV.

    Args:
    ----
        epoch (Union[int, None])             : Current epoch.
        lr (Union[float, None])              : Learning rate.
        training_loss (Union[float, None])   : Training loss.
        validation_loss (Union[float, None]) : Validation loss.
        score (Union[float, None])           : Model score (the UC-JSD).
        append (bool)                        : If True, appends to existing file.
                                               If False, creates a new file.
    """
    convergence_path = constants.job_dir + "convergence.log"
    if constants.job_type == "fine-tune":
        epoch_label = "Step"
    else:
        epoch_label = "Epoch"

    if not append:  # create the file
        with open(convergence_path, "w") as output_file:
            # write the header
            output_file.write(f"{epoch_label.lower()}, lr, avg_train_loss, "
                              f"avg_valid_loss, model_score\n")
    else:  # append to existing file
        # only write a `convergence.log` when training
        if constants.job_type in ["train", "fine-tune"]:
            if score is None:
                with open(convergence_path, "a") as output_file:
                    output_file.write(f"{epoch_label} {epoch}, {lr:.8f}, "
                                      f"{training_loss:.8f}, "
                                      f"{validation_loss:.8f}, ")
                # write to tensorboard
                tb_writer.add_scalar("Training/training_loss", training_loss, epoch)
                tb_writer.add_scalar("Training/validation_loss", validation_loss, epoch)
                tb_writer.add_scalar("Training/lr", lr, epoch)

            elif score == "NA":
                with open(convergence_path, "a") as output_file:
                    output_file.write(f"{score}\n")

            elif score is not None and training_loss is not None:
                with open(convergence_path, "a") as output_file:
                    output_file.write(f"{epoch_label} {epoch}, {lr:.8f}, "
                                      f"{training_loss:.8f}, "
                                      f"{validation_loss:.8f}, {score:.6f}\n")

            elif score is not None:
                with open(convergence_path, "a") as output_file:
                    output_file.write(f"{score:.6f}\n")

            else:
                raise NotImplementedError

def write_molecules(molecules : list,
                    final_likelihoods : torch.Tensor,
                    epoch : str,
                    write : bool=False,
                    label : str="test") -> Tuple[list, list, list]:
    """
    Writes generated molecular graphs and their NLLs. In writing the structures
    to a SMILES file, determines if structures are valid and returns this
    information (to avoid recomputing later).

    Args:
    ----
        molecules (list)                 : Contains generated `MolecularGraph`s.
        final_likelihoods (torch.Tensor) : Contains final NLLs for the graphs.
        epoch (str)                      : Number corresponding to the current
                                           training epoch.
        write (bool)                     : Whether or not to write the molecules
                                           to SMILES.
        label (str)                      : Label to use when saving the molecules.

    Returns:
    -------
        fraction_valid (float)           : The fraction of molecular graphs which are
                                           valid molecules.
        validity_tensor (torch.Tensor)   : A binary vector indicating the chemical
                                           validity of the input graphs using either
                                           a 1 (valid) or 0 (invalid).
        uniqueness_tensor (torch.Tensor) : A binary vector indicating the uniqueness
                                           of the input graphs using either a 1
                                           (unique or first duplicate instance) or
                                           0 (duplicate).
    """
    # save molecules as SMILE
    if constants.job_type == "fine-tune":
        step         = epoch.split(" ")[1]
        smi_filename = constants.job_dir + f"generation/step{step}_{label}.smi"
    else:
        epoch_number = epoch.split(" ")[1]
        smi_filename = constants.job_dir + f"generation/epoch_{epoch_number}.smi"

    (fraction_valid,
     validity_tensor,
     uniqueness_tensor) = write_graphs_to_smi(smi_filename=smi_filename,
                                              molecular_graphs_list=molecules,
                                              write=write)
    # save the NLLs and validity status
    write_likelihoods(likelihood_filename=f"{smi_filename[:-3]}likelihood",
                      likelihoods=final_likelihoods)
    write_validity(validity_file_path=f"{smi_filename[:-3]}valid",
                   validity_tensor=validity_tensor)

    return fraction_valid, validity_tensor, uniqueness_tensor

def write_likelihoods(likelihood_filename : str,
                      likelihoods : torch.Tensor) -> None:
    """
    Writes the final likelihoods of each molecule to a file in the same order as
    the molecules are written in the corresponding SMILES file.

    Args:
    ----
        likelihood_filename (str)  : Path to the likelihood-containing file.
        likelihoods (torch.Tensor) : Likelihoods.
    """
    with open(likelihood_filename, "w") as likelihood_file:
        for likelihood in likelihoods:
            likelihood_file.write(f"{likelihood}\n")

def write_ts_properties(training_set_properties : dict) -> None:
    """
    Writes the training set properties to CSV.

    Args:
    ----
        training_set_properties (dict) : The properties of the training set.
    """
    training_set = constants.training_set  # path to "train.smi"
    dict_path    = f"{training_set[:-4]}.csv"

    with open(dict_path, "w") as csv_file:

        csv_writer = csv.writer(csv_file, delimiter=";")
        for key, value in training_set_properties.items():
            if "validity_tensor" in key:
                # skip writing the validity tensor here because it is really
                # long, instead it gets its own file elsewhere
                continue
            if isinstance(value, np.ndarray):
                csv_writer.writerow([key, list(value)])
            elif isinstance(value, torch.Tensor):
                try:
                    csv_writer.writerow([key, float(value)])
                except ValueError:
                    csv_writer.writerow([key, [float(i) for i in value]])
            else:
                csv_writer.writerow([key, value])

def write_validation_scores(output_dir : str, epoch_key : str,
                            model_scores : dict, append : bool=True) -> None:
    """
    Writes a CSV with the model validation scores as a function of the epoch.

    Args:
    ----
        output_dir (str)    : Full path/filename to CSV file.
        epoch_key (str)     : For example, "Training set" or "Epoch {n}".
        model_scores (dict) : Contains the average NLL per molecule of
                              {validation/train/generated} structures, and the
                              average model score (weighted mean of above two scores).
        append (bool)       : Indicates whether to append to the output file or
                              start a new one.
    """
    validation_file_path = output_dir + "validation.log"
    avg_likelihood_val   = model_scores["avg_likelihood_val"]
    avg_likelihood_train = model_scores["avg_likelihood_train"]
    avg_likelihood_gen   = model_scores["avg_likelihood_gen"]
    uc_jsd               = model_scores["UC-JSD"]

    if not append:  # create file
        with open(validation_file_path, "w") as output_file:
            # write headeres
            output_file.write("set, avg_likelihood_per_molecule_val, "
                              "avg_likelihood_per_molecule_train, "
                              "avg_likelihood_per_molecule_gen, uc_jsd\n")

    # append the properties of interest to the CSV file
    with open(validation_file_path, "a") as output_file:
        output_file.write(f"{epoch_key:}, {avg_likelihood_val:.5f}, "
                          f"{avg_likelihood_train:.5f}, "
                          f"{avg_likelihood_gen:.5f}, {uc_jsd:.7f}\n")

    try:  # write to tensorboard
        epoch = int(epoch_key.split()[1])
        # scalars
        tb_writer.add_scalar("Evaluation/avg_validation_likelihood", avg_likelihood_val, epoch)
        tb_writer.add_scalar("Evaluation/avg_training_likelihood", avg_likelihood_train, epoch)
        tb_writer.add_scalar("Evaluation/avg_generation_likelihood", avg_likelihood_gen, epoch)
        tb_writer.add_scalar("Evaluation/uc_jsd", uc_jsd, epoch)
    except:
        pass

def write_validity(validity_file_path : str,
                   validity_tensor : torch.Tensor) -> None:
    """
    Writes the validity (0 or 1) of each molecule to a file in the same
    order as the molecules are written in the corresponding SMILES file.

    Args:
    ----
        validity_file_path (str)       : Path to validity file.
        validity_tensor (torch.Tensor) : A binary vector indicating the chemical
                                         validity of the input graphs using either
                                         a 1 (valid) or 0 (invalid).
    """
    with open(validity_file_path, "w") as valid_file:
        for valid in validity_tensor:
            valid_file.write(f"{valid}\n")

def tbwrite_loglikelihoods(step : Union[int, None]=None,
                           agent_loglikelihoods : Union[torch.Tensor, None]=None,
                           prior_loglikelihoods : Union[torch.Tensor, None]=None) -> \
                           None:
    """
    Writes the current epoch and log-likelihoods to the tensorboard during
    fine-tuning jobs.

    Args:
    ----
        step (Union[int, None])                          : Current fine-tuning step.
        agent_loglikelihoods (Union[torch.Tensor, None]) : Vector containing the
                                                           agent log-likelihoods
                                                           for sampled molecules.
        prior_loglikelihoods (Union[torch.Tensor, None]) : Vector containing the
                                                           prior log-likelihoods
                                                           for sampled molecules.
    """
    avg_agent_loglikelihood = torch.mean(agent_loglikelihoods)
    avg_prior_loglikelihood = torch.mean(prior_loglikelihoods)
    tb_writer.add_scalar("Train/agent_loglikelihood", avg_agent_loglikelihood, step)
    tb_writer.add_scalar("Train/prior_loglikelihood", avg_prior_loglikelihood, step)

def load_saved_model(model : torch.nn.Module, path : str) -> torch.nn.Module:
    """
    Loads a pre-saved neural net model.

    Args:
    ----
        model (torch.nn.Module) : Existing but bare-bones model variable.
        path (str)              : Path to the saved model.

    Returns:
    -------
        model (torch.nn.Module) : Loaded model.
    """
    try:
        # first try to load model as if it was created using GraphINVENT
        # v1.0 (will raise an exception if it was actually created with
        # GraphINVENT v2.0)
        model.state_dict = torch.load(path).state_dict()
    except AttributeError:
        # try to load the model as if created using GraphINVENT v2.0 or
        # later
        model.load_state_dict(torch.load(path))
    return model

max_stars_repo_name: nicoBruhn/GraphINVENT

max_issues_repo_name: nicoBruhn/GraphINVENT

max_forks_repo_name: nicoBruhn/GraphINVENT


text:
# reimplementation of integral_2d1d_v01 for up to date bundles

from load import ROOT as R
import numpy as N
import gna.constructors as C
from gna.bundle import TransformationBundle

class integral_2d1d_v02(TransformationBundle):
    def __init__(self, *args, **kwargs):
        TransformationBundle.__init__( self, *args, **kwargs )
        self.check_cfg()

    @staticmethod
    def _provides(cfg):
        var0, var1 = cfg.variables
        return (), ('integral', var0, var1, var0+'_edges', var0+'_hist')

    def check_cfg(self):
        if not 'name' in self.cfg:
            pkey = self.cfg.parent_key()
            if not pkey:
                raise Exception('"name" option is not provided for integral_1d_v01')
            self.cfg.name = pkey

        try:
            self.edges = N.ascontiguousarray(self.cfg.edges, dtype='d')
        except:
            raise Exception('Invalid binning definition: {!r}'.format(self.cfg.edges))

        try:
            self.xorders = N.ascontiguousarray(self.cfg.xorders, dtype='P')
        except:
            raise Exception('Invalid xorders definition: {!r}'.format(self.cfg.xorders))

        if len(self.cfg.variables)!=2:
            raise Exception('Two vairables should be provided')

    def build(self):
        if self.xorders.size>1:
            if self.xorders.size+1 != self.edges.size:
                raise Exception('Incompartible edges and xorders definition:\n    {!r}\n    {!r}'.format(self.edges, self.xorders))
            self.integrator = R.GaussLegendre2d(self.edges, self.xorders, self.edges.size-1, -1.0, 1.0, self.cfg.yorder)
        else:
            self.integrator = R.GaussLegendre2d(self.edges, int(self.xorders[0]), self.edges.size-1, -1.0, 1.0, self.cfg.yorder)
        self.integrator.points.setLabel('Gauss-Legendre 2d')

        self.integrator.points.x.setLabel(self.cfg.variables[0])
        self.integrator.points.xedges.setLabel('%s edges'%self.cfg.variables[0])
        self.integrator.points.y.setLabel(self.cfg.variables[1])

        self.set_output(self.cfg.variables[0],            None, self.integrator.points.x)
        self.set_output('%s_edges'%self.cfg.variables[0], None, self.integrator.points.xedges)
        self.set_output('%s_hist'%self.cfg.variables[0],  None, self.integrator.points.xhist)
        self.set_output(self.cfg.variables[1],            None, self.integrator.points.y)

        for i, it in enumerate(self.nidx):
            hist = R.GaussLegendre2dHist(self.integrator)
            hist.hist.setLabel(it.current_format(name='hist'))

            self.set_input('integral', it, hist.hist.f, argument_number=0)
            self.set_output('integral', it, hist.hist.hist)

    def define_variables(self):
        pass

ext: py

max_stars_repo_name: gnafit/gna

max_issues_repo_name: gnafit/gna

max_forks_repo_name: gnafit/gna

reason: import numpy


text:
"""Poisson equation."""

import numpy as np
from skfem.assembly import bilinear_form, linear_form

@bilinear_form
def laplace(u, du, v, dv, w):
    return sum(du * dv)

@bilinear_form
def vector_laplace(u, du, v, dv, w):
    from .helpers import ddot
    return ddot(du, dv)

@bilinear_form
def mass(u, du, v, dv, w):
    return u * v

@linear_form
def unit_load(v, dv, w):
    return v


text:
############################################################################################
# Program Name : Confidence Interval Calculator / Grapher
# Description  : Displays confidence intervals from actual MOS values and calculated
#                   statistics.  Made from a matlab version by a video quality researcher
#                   in germany.
############################################################################################

import os
import sys
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import scipy
import scipy.io

# Input Value Files and Fields
MOSFileNameList     = []
MOSFieldNameList    = []
NRParsFileNameList  = []
NRParsFieldNameList = []
NRParsMOSCount      =  0

GraphSaveFileName =    ""
IsVerbose         = False

# Global field values
MetricName = ""
MOSDict    = []
NRParsDict = []

# Main
#   Plot the confidence interval (CI) of an NR parameter by MOS values
# SYNTAX
#   python ci_calc.py -m iqa_camera.mat ccriq_dataset blur.mat 10_percent
# SEMANTICS
#   Read values and run ci_calc algorithm. Please see ci_calc
#   description or documentation for extra detail on f output graph
#   and output parameters.
#
# Input Parameters:
#   mosFileName     Input mos filename(s)
#   mosFieldName    MOS filename's field name
#   nrParsFileName  Input NR Parameter's filename(s)
#   nrParsFieldName NR Parameter's field name
##
# Output Parameters
#   ideal_ci = the ideal confidence interval
#   practial_ci = the practical confidence interval
#
# Constraints:
#   All datasets are weighted equally.
#   The MOSs must range from 1 to 5.
#
def Main():
    # Parse commandline arguments and read mat files
    parse_command_arguments(sys.argv)
    read_mos_and_nrpars()

    # Calculate confidence intervals
    (ideal_ci, practical_ci) = ci_calc(MetricName, MOSDict, NRParsDict, GraphSaveFileName, IsVerbose)

# Main
#   Parse command arguments
# SYNTAX
#   parse_command_arguments([ "-help" ])
# SEMANTICS
#   Parses input filenames and fields.  And, prints any format errors.
#
# Input Parameters:
#   argvList: Command line arguments
#
# Output Parameters
#   MOSFileNameList     Input MOS filename(s)
#   MOSFieldNameList    MOS filename's field name
#   NRParsFileNameList  Input NR Parameter's filename(s)
#   NRParsFieldNameList NR Parameter's field name
#
def parse_command_arguments(argvList):
    # Global vars
    global MOSFileNameList
    global MOSFieldNameList
    global NRParsFileNameList
    global NRParsFieldNameList
    global NRParsMOSCount
    global GraphSaveFileName
    global IsVerbose

    # Clear data
    MOSFileNameList     =    []
    MOSFieldNameList    =    []
    NRParsFileNameList  =    []
    NRParsFieldNameList =    []
    NRParsMOSCount      =     0
    GraphSaveFileName   =    ""
    IsVerbose           = False

    # Parse arguments
    if len(argvList)<=1:
        argvList = [ argvList[0] if len(argvList)>0 else "", "-h" ]

    argvListCount = len(argvList)
    index = 1

    while index<argvListCount:
        # Format argument
        argv = (argvList[index] or "").strip().lower()
        if argv.startswith("/"):
            argv = "-"+argv[1:]

        # Read argument
        if (argv=="-m" and index+4<argvListCount):
            if NRParsMOSCount>0:
                print("  <Warning: Only 1 Dataset Can Currently be Plotted, Defaulting to Last Entered Dataset>")
            MOSFileNameList    .append((argvList[index+1] or "").strip())
            MOSFieldNameList   .append((argvList[index+2] or "").strip())
            NRParsFileNameList .append((argvList[index+3] or "").strip())
            NRParsFieldNameList.append((argvList[index+4] or "").strip())
            NRParsMOSCount += 1
            index += 4
        elif (argv=="-s" and index+1<argvListCount):
            GraphSaveFileName = (argvList[index+1] or "").strip()
            index += 1
        elif (argv=="-h" or argv=="--help"):
            print("Usage:")
            print(" python3 ci_calc.py [options]")
            print("Options:")
            print(" -m    <mosFileName> <mosFieldName> <nrParsFileName> <nrParsFieldName>")
            print(" -s    <graphSaveFileName>")
            print(" -b    Verbose / Program Status")
            print("Misc Options:")
            print(" -h --help       Help")
            print(" -v --version    Version Number")
            print(" -b --verbose    Verbose Messages")
            print("Example:")
            print(" python3 ci_calc.py -b -m iqa_camera.mat ccriq_dataset \"..\\nr_data\\group_blur\\NRpars_blur_ccriq.mat\" unsharp -s \"ci_graph1.jpg\"")
            print(" python3 ci_calc.py -b -m iqa_camera.mat ccriq_dataset \"..\\nr_data\\group_blur\\NRpars_blur_ccriq.mat\" viqet-sharpness -s \"ci_graph1.jpg\"")
            print(" python3 ci_calc.py -v")
        elif (argv=="-v" or argv=="--version"):
            print("Version 1.0a")
        elif (argv=="-b" or argv=="-?" or argv=="--verbose"):
            IsVerbose = True
        else:
            print("  <Error: Failed to Parse Arguments at Index {0} of {1} with First Value \"{2}\">".format(index, argvListCount, argvList[index]))

        index += 1

# read_mos_and_nrpars
#   Read MOS and NRPars
# SYNTAX
#   read_mos_and_nrpars();
# SEMANTICS
#   Read MOS and NRPars into global variables.
#
# Input Parameters:
#   MOSFileNameList       MOS filename(is a list)
#   MOSFieldNameList      MOS dataset name
#   NRParsFileNameList    NRPars filename
#   NRParsFieldNameList   NRPars fieldname
#
# Output Parameters
#   MetricName   Metric name, taken from fieldnames
#   MOSDict      MOS values from mat file
#   NRParsDict   NRPar values from mat file
#
def read_mos_and_nrpars():
    # Global vars
    global MetricName
    global MOSDict
    global NRParsDict

    MetricName = ""
    MOSDict    = {}
    NRParsDict = {}

    # Print header
    if IsVerbose:
        print("Reading MOS and NRPars from Files")

    # Set metric name
    for index in range(NRParsMOSCount):
        MetricName += "{0}({1}) ".format(MOSFieldNameList[index], NRParsFieldNameList[index])
    MetricName = MetricName.strip()

    # Read mat files
    for index in range(NRParsMOSCount):
        # Print header
        if IsVerbose:
            print("Reading Files \"{0}\" and \"{1}\"".format(MOSFileNameList[index], NRParsFileNameList[index]))

        # Set MOS and NRPar names
        dataSetName = "{0}({1}) ".format(MOSFileNameList[index], MOSFieldNameList[index], NRParsFileNameList[index], NRParsFieldNameList[index])

        # Attempt to read MOS file
        if not os.path.exists(MOSFileNameList[index]):
            print("  <Failed to find MOS File \"{0}\">".format(MOSFileNameList[index]))
            exit(0)

        data = scipy.io.loadmat(MOSFileNameList[index])

        if not MOSFieldNameList[index] in data.keys():
            print("  <Failed to find MOS Dataset Field \"{0}\">".format(MOSFieldNameList[index]))
            exit(0)

        data = data[MOSFieldNameList[index]]

        if not( type(data      )==np.ndarray and len(data   )>0 \
            and type(data[0]   )==np.ndarray and len(data[0])>0 \
            and type(data[0][0])==np.void    and "media" in np.dtype(data[0][0]).names):
            print("  <Failed to find \"media\" Field in MOS Dataset \"{0}\">".format(MOSFieldNameList[index]))
            exit(0)

        dataIndex = np.dtype(data[0][0]).names.index("media")
        data = data[0][0]

        if not(dataIndex<len(data) \
            and type(data[dataIndex]   )==np.ndarray and len(data[dataIndex]   )>0 \
            and type(data[dataIndex][0])==np.ndarray and len(data[dataIndex][0])>0):
            print("  <Failed to find \"media\" Field as a List in MOS Dataset \"{0}\">".format(MOSFieldNameList[index]))
            exit(0)

        MOSDict[dataSetName] = []
        for rowIndex in range(len(data[dataIndex][0])):
            row = data[dataIndex][0][rowIndex]
            if not("mos" in np.dtype(row).names):
                print("  <Failed to find MOS column on row \"{0}\">".format(rowIndex))
                exit(0)
            columnIndex = np.dtype(row).names.index("mos")
            if not( type(row[columnIndex]      )==np.ndarray and len(row[columnIndex]   )>0 \
                and type(row[columnIndex][0]   )==np.ndarray and len(row[columnIndex][0])>0 \
                and type(row[columnIndex][0][0]) in (np.int8, np.int16, np.int32, np.int64, np.uint8, np.uint16, np.uint32, np.uint64, np.float32, np.float64)):
                print("  <Failed to find MOS column as a formatted float on row \"{0}\">".format(rowIndex))
                exit(0)
            mos = float(row[columnIndex][0][0])
            MOSDict[dataSetName].append(mos)

        # Attempt to read NRPars file
        if not os.path.exists(NRParsFileNameList[index]):
            print("  <Failed to find NRPars File \"{0}\">".format(NRParsFileNameList[index]))
            exit(0)

        data = scipy.io.loadmat(NRParsFileNameList[index])

        if not "NRpars" in data.keys():
            print("  <Failed to find NRPars Dataset Field \"NRPars\">")
            exit(0)

        data = data["NRpars"]

        if not( type(data      )==np.ndarray and len(data      )>0 \
            and type(data[0]   )==np.ndarray and len(data[0]   )>0 \
            and type(data[0][0])==np.void    and "par_name" in np.dtype(data[0][0]).names and "data" in np.dtype(data[0][0]).names):
            print("  <Failed to find \"par_name\" and \"data\" Fields in NRPars Dataset \"{0}\">".format(NRParsFieldNameList[index]))
            exit(0)

        data = data[0][0]
        dataIndex1 = np.dtype(data).names.index("par_name")
        dataIndex2 = np.dtype(data).names.index("data")

        if not(dataIndex1<len(data) and dataIndex2<len(data) \
            and type(data[dataIndex1]   )==np.ndarray and len(data[dataIndex1]   )>0 \
            and type(data[dataIndex1][0])==np.ndarray and len(data[dataIndex1][0])>0 \
            and type(data[dataIndex2]   )==np.ndarray and len(data[dataIndex2]   )>0):
            print("  <Failed to find \"par_name\" and \"data\" Fields as a List in MOS Dataset \"{0}\">".format(MOSFieldNameList[index]))
            exit(0)

        parNames = [ str(d[0]) for d in data[dataIndex1][0] if type(d)==np.ndarray and len(d)>0 and type(d[0])==np.str_ ]

        if not(NRParsFieldNameList[index] in parNames):
            print("  <Failed to find parameter name \"{0}\">".format(NRParsFieldNameList[index]))
            exit(0)

        parNamesIndex = parNames.index(NRParsFieldNameList[index])

        if not(parNamesIndex<len(data[dataIndex2]) \
            and type(data[dataIndex2][parNamesIndex])==np.ndarray and len(data[dataIndex2][parNamesIndex])>0):
            print("  <Failed to find data Fields for parameter name \"{0}\">".format(NRParsFieldNameList[index]))
            exit(0)

        NRParsDict[dataSetName] = []
        for rowIndex in range(len(data[dataIndex2][parNamesIndex])):
            val = data[dataIndex2][parNamesIndex][rowIndex]
            if not type(val) in (np.int8, np.int16, np.int32, np.int64, np.uint8, np.uint16, np.uint32, np.uint64, np.float32, np.float64):
                print("  <Failed to find data as a formatted float on row \"{0}\">".format(rowIndex))
                exit(0)
            NRParsDict[dataSetName].append(val)

        if not len(MOSDict[dataSetName])==len(NRParsDict[dataSetName]):
            print("  <Failed To Read Data Where {0} MOS Rows and {1} NRPars Rows>".format(len(MOSDict[dataSetName]), len(NRParsDict[dataSetName])))
            exit(0)

    # Print footer
    if IsVerbose:
        print("Reading MOS and NRPars from Files Complete")

# ci_calc
#   Estimate the confidence interval (CI) of an NR parameter
# SYNTAX
#   (ideal_ci, practical_ci) = ci_calc(metric_name, dataset_mos,
#       dataset_metrics, fig_path = False, verbose = True);
# SEMANTICS
#   Estimate the confidence interval (CI) of an NR metric or parameter,
#   by comparing the conclusions reached by the metric with conclusions
#   reached by a subjective test. Both will use a constant confidence
#   interval (CI) to make decisions. The subjective CI is based on
#   5-level ACR MOSs. Two recommended CIs are printed to the command window.
#   (1) ideal CI, and (2) practical CI. The classification types are plotted,
#   which allows the user to choose an alternate CI.
#
#   By analogy, assess the performance of the metric in terms of an ad-hoc
#   test with N people. This analysis assumes that the metric and MOSs are
#   compared without statistical tests or confidence intervals.
#
# Input Parameters:
#   metric_name     Character string that contains the metric's name
#   dataset_mos     Cell array. For each dataset (1..num_datasets), a
#                   double array that contains the mean opinion score (MOS)
#                   for each stimuli in the dataset.
#   dataset_metrics Cell array. For each dataset (1..num_datasets), a
#                   double array that contains the metric's value for each
#                   stimuli in the dataset. Order of stimuli must be
#                   identical to dataset_mos.
#   fig_path        figure path for saving
#   verbose         print extra status messages
#
#   The theoretical underpinnings of this algorithm are pending publication
#   of NTIA Report "Confidence Intervals for Subjective Tests and
#   Objective Metrics" by Margaret H Pinson
#
#   For a preliminary analysis, see Margaret Pinson, "NR metric confidence
#   interval estimation using classification errors," Video Quality Experts
#   Group (VQEG) meeting, Statistical Analysis Methods (SAM) Group,
#   Presentation 11, March 2020.
#   ftp://vqeg.its.bldrdoc.gov/Documents/VQEG_online_Mar20/VQEG_2020_SAM_011_confidence_intervals_for_metrics.pptx
#
# Output Parameters
#   ideal_ci = the ideal confidence interval
#   practial_ci = the practical confidence interval
#   N = the number of people in an ad-hoc test with an equivalent likelihood of
#       false ranking, or zero (0) if the performance is worse than a 1
#       person ad-hoc test.
#
#   For positively correlated metrics, false ranking is where a well designed
#   subjective test would conclude that stimuli A is statistical better
#   than stimuli B, but the metric value for stimuli B is greater than
#   the metric value for stimuli A. "Less than" is used for negatively
#   correlated metrics.
#
# Constraints:
#   All datasets are weighted equally.
#   The MOSs must range from 1 to 5.
#
def ci_calc(metric_name, dataset_mos, dataset_metrics, fig_path = False, verbose = True):
    threshold_level = 0.5 # delta S, where 95% of stimuli MOS can be rank ordered
    false_rank_thresh = 0.01 # disagree rate
    false_diff_thresh = 0.10 # half of the uncertain rate of 20%
    practical_threshold = 0.165 # half of maximum uncertain rate plus disagree rate
    concur_threshold = 0.91 #based on analyses of the VQEG FRTV Phase I ratings
    print('Metric confidence interval analysis for {}'.format(metric_name))

    # calculate range of this parameter
    this_par = []
    pos_corr = []
    for dcnt, mos in dataset_mos.items():
        omos = dataset_metrics[dcnt]
        tmp = np.corrcoef(mos, omos)
        if tmp[0,1] >= 0:
            pos_corr.append(1)
        else:
            pos_corr.append(-1)
        this_par.extend(omos)

    this_par = sorted(this_par)
    this_par_count = len(this_par)
    pmin = min(this_par)
    pmax = max(this_par)

    if verbose and len(this_par)>0:
        print('Full range {}..{}, '.format(pmin, pmax))
        print('95% of data in {}..{}\n'.format(this_par[max(int(round(0.025*this_par_count)), this_par_count-1)], this_par[max(int(round(0.975*this_par_count)), this_par_count-1)]))

    if sum(pos_corr) > 0:
        if verbose:
            print('Positively correlated with MOS for most datasets\n\n')
        is_pos_corr = True
    elif sum(pos_corr) == 0:
        if verbose:
            print('Split decision on whether metric is positively or negatively correlated with MOS.\nAssume positive correlation.\n\n')
        is_pos_corr = True
    else:
        if verbose:
            print('Negatively correlated with MOS for most datasets\n\n')
        is_pos_corr = False

    if pmin == pmax:
        print('Warning: parameter has a constant value, aborting.\n')
        ideal_ci = np.NaN
        practical_ci = np.NaN
        return

    # manually loop through all pairs of stimuli
    subj = []
    obj = []
    wt = []
    curr = 1
    for dcnt, mos in dataset_mos.items():
        curr_len = len(mos)
        for mcnt1 in range(0,curr_len):
            for mcnt2 in range(mcnt1+1,curr_len):
                # subj(curr) is decision whether #1 is better,
                # equivalent, or worse than #2
                diff = mos[mcnt1] - mos[mcnt2]
                if diff > threshold_level:
                    subj.append(1)
                elif diff < -threshold_level:
                    subj.append(-1)
                else:
                    subj.append(0)

                # obj(curr) is distance before thresholding, since the
                # point of this function is to ideal_ci a threshold
                obj.append(dataset_metrics[dcnt][mcnt1] - dataset_metrics[dcnt][mcnt2])

                # note weight
                wt.append( 1 / len(mos) )
                curr += 1

    # flip sign of objective differences, if parameter is
    # negatively correlated to MOS
    if is_pos_corr == False:
        obj = list(-np.array(obj))

    # Have all of the data. Now make the plot.
    # round our increment to one significant digits
    # incr = round((pmax-pmin)/100, 1, 'significant');
    incr = (pmax-pmin)/100
    list_want = np.arange(incr,pmax-pmin,incr)

    correct_rank = [0]*len(list_want)
    correct_tie = [0]*len(list_want)
    false_ranking = [0]*len(list_want)
    false_distinction = [0]*len(list_want)
    false_tie = [0]*len(list_want)

    # create data for roughly 60% of the range of parameter values
    # from there, the plot flattens and contains no more info
    for loop in range(0,len(list_want)):
        #print(loop)
        delta = list_want[loop]
        for curr in range(0,len(subj)):
            if (subj[curr] == 1 and obj[curr] >= delta) or (subj[curr] == -1 and obj[curr] <= -delta):
                correct_rank[loop] = correct_rank[loop] + wt[curr]
            elif subj[curr] == 0 and obj[curr] > -delta and obj[curr] < delta:
                correct_tie[loop] = correct_tie[loop] + wt[curr]
            elif (subj[curr] == 1 and obj[curr] <= -delta) or (subj[curr] == -1 and obj[curr] >= delta):
                false_ranking[loop] = false_ranking[loop] + wt[curr]
            elif (subj[curr] != 0 and obj[curr] > -delta and obj[curr] < delta):
                false_tie[loop] = false_tie[loop] + wt[curr]
            else:
                false_distinction[loop] = false_distinction[loop] + wt[curr]

    total_votes = sum(wt)

    correct_rank = [ val/total_votes for val in correct_rank]
    correct_tie = [ val/total_votes for val in correct_tie]
    false_ranking = [ val/total_votes for val in false_ranking]
    false_distinction = [ val/total_votes for val in false_distinction]
    false_tie = [ val/total_votes for val in false_tie]

    # if too much data is false_tie and correct_tie at minimum
    # threshold, don't try. Skip. Rule of thumb: 50% ties. We expect
    # values close to zero, so this should mean most of the metric is a
    # constant value.
    if false_tie[0] + correct_tie[0] > 0.5:
        print('Half of data is correct ties or false ties. Skipping.\n')
        ideal_ci = np.NAN
        practical_ci = np.NAN
        return

    # compute the ideal ci
    ideal_ci = len(list_want)
    for n, (fr, fd) in enumerate(zip(false_ranking,false_distinction)):
        if (fr < false_rank_thresh) and (fd < false_diff_thresh):
            ideal_ci = n
            break

    # compute the practical CI
    practical_ci = len(list_want)
    for n, (fr, fd) in enumerate(zip(false_ranking, false_distinction)):
        if (fr + fd) < practical_threshold:
            practical_ci = n
            break

    equiv_ideal = np.sqrt(correct_rank[ideal_ci]) + 1.2 * correct_tie[ideal_ci]
    equiv_practical = np.sqrt(correct_rank[practical_ci]) + 1.2 * correct_tie[practical_ci]

    # print recommended threshold
    if verbose:
        print('{} Ideal CI      ({}% correct ranking, {}% false ranking, {}% false distinction, {}% false tie, {}% correct tie)'.format( \
            list_want[ideal_ci], round(correct_rank[ideal_ci]*100), round(false_ranking[ideal_ci]*100), round(false_distinction[ideal_ci]*100), \
            round(false_tie[ideal_ci]*100), round(correct_tie[ideal_ci]*100)))
        if equiv_ideal >= concur_threshold:
            print(' ==> equivalent to a subjective test with 24 subjects')

        print('\n{} Practical CI  ({}% correct ranking, {}% false ranking, {}% false distinction, {}% false tie, {}% correct tie)'.format( \
            list_want[practical_ci], round(correct_rank[practical_ci]*100), round(false_ranking[practical_ci]*100), round(false_distinction[practical_ci]*100), \
            round(false_tie[practical_ci]*100), round(correct_tie[practical_ci]*100)))
        if equiv_practical >= concur_threshold:
            print(' ==> equivalent to a subjective test with 15 subjects')

    # create plot
    fig, ax = plt.subplots(figsize=(6, 6))
    ax.plot(list_want, [val*100 for val in correct_rank], 'g', label="correct rank")
    ax.plot(list_want, [val*100 for val in false_ranking], 'r', label="false rank")
    ax.plot(list_want, [val*100 for val in false_distinction], 'b--', label="false distinction")
    ax.plot(list_want, [val*100 for val in false_tie], 'y--', label="false tie")
    ax.plot(list_want, [val*100 for val in correct_tie], 'y', label="correct tie")

    ax.plot([list_want[ideal_ci],list_want[ideal_ci]], [0,100], 'k', label="ideal CI")
    ax.plot([list_want[practical_ci],list_want[practical_ci]], [0,100], 'k--', label="practical CI")

    plt.xlabel("Delta Metric")
    plt.ylabel("Probability")
    ax.grid()
    ax.set_title(metric_name)
    ax.axis([0, 2, 0, 100])
    leg = ax.legend()

    if fig_path:
        plt.savefig(fig_path)
    else:
        plt.show()

    # equivelence determination

    correct_rank_zero = 0
    correct_tie_zero = 0
    false_ranking_zero = 0
    false_distinction_zero = 0
    false_tie_zero = 0

    delta = 0
    for curr in range(0,len(subj)):
        if (subj[curr] == 1 and obj[curr] >= delta) or (subj[curr] == -1 and obj[curr] <= -delta):
            correct_rank_zero = correct_rank_zero + wt[curr]
        elif subj[curr] == 0 and obj[curr] > -delta and obj[curr] < delta:
            correct_tie_zero = correct_tie_zero + wt[curr]
        elif (subj[curr] == 1 and obj[curr] <= -delta) or (subj[curr] == -1 and obj[curr] >= delta):
            false_ranking_zero = false_ranking_zero + wt[curr]
        elif (subj[curr] != 0 and obj[curr] > -delta and obj[curr] < delta):
            false_tie_zero = false_tie_zero + wt[curr]
        else:
            false_distinction_zero = false_distinction_zero + wt[curr]

    correct_rank_zero = correct_rank_zero/total_votes
    correct_tie_zero = correct_tie_zero/total_votes
    false_ranking_zero = false_ranking_zero/total_votes
    false_distinction_zero = false_distinction_zero/total_votes
    false_tie_zero = false_tie_zero/total_votes

    print('\nNo CI used ({}% correct ranking, {}% false ranking, {}% false distinction, {}% false tie, {}% correct tie)'.format( \
        round(correct_rank_zero*100), round(false_ranking_zero*100), round(false_distinction_zero*100), round(false_tie_zero*100), round(correct_tie_zero*100)))

    if false_ranking_zero <= 0.0325:
        equivalent = 12
        print(' ==> equivalent to a pilot test with {} subjects'.format(equivalent))
    elif false_ranking_zero <= 0.0395:
        equivalent = 9
        print(' ==> equivalent to a pilot test with {} subjects'.format(equivalent))
    elif false_ranking_zero <= 0.056:
        equivalent = 6
        print(' ==> equivalent to a pilot test with {} subjects'.format(equivalent))
    elif false_ranking_zero <= 0.0765:
        equivalent = 3
        print(' ==> equivalent to a {} person ad-hoc test'.format(equivalent))
    elif false_ranking_zero <= 0.0995:
        equivalent = 2
        print(' ==> equivalent to a {} person ad-hoc test'.format(equivalent))
    elif false_ranking_zero <= 0.1285:
        equivalent = 1
        print(' ==> equivalent to a {} person ad-hoc test'.format(equivalent))
    else:
        equivalent = 0

    return list_want[ideal_ci], list_want[practical_ci]

if __name__ == "__main__":
    Main()

max_stars_repo_name: NTIA/NRMetricFramework

max_issues_repo_name: NTIA/NRMetricFramework

max_forks_repo_name: NTIA/NRMetricFramework

reason: import numpy,import scipy


text:
# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2015-2018 GEM Foundation
#
# OpenQuake is free software: you can redistribute it and/or modify it
# under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# OpenQuake is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>.

"""
Module exports :class:`DrouetBrazil2015`
               :class:`DrouetBrazil2015_with_depth`
"""
import numpy as np

from openquake.hazardlib.gsim.base import CoeffsTable, GMPE
from openquake.hazardlib import const
from openquake.hazardlib.imt import PGA, PGV, SA
from scipy.constants import g

class DrouetBrazil2015(GMPE):
    """
    Implements GMPE developed by S. Drouet unpublished for Brazil based on the
    method described in Douet & Cotton (2015) BSSA doi: 10.1785/0120140240.
    """

    #: Supported tectonic region type is stable continental crust given that
    #: the equations have been derived for Eastern North America.
    DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.STABLE_CONTINENTAL

    #: Supported intensity measure types are spectral acceleration,
    #: and peak ground acceleration, see table 6, page 1022 (PGA is assumed
    #: to be equal to SA at 0.01 s)
    DEFINED_FOR_INTENSITY_MEASURE_TYPES = set([
        PGA,
        PGV,
        SA
    ])

    #: Supported intensity measure component is the geometric mean of
    #two : horizontal components
    #:attr:`~openquake.hazardlib.const.IMC.AVERAGE_HORIZONTAL`,
    DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.AVERAGE_HORIZONTAL

    #: Supported standard deviation type is only total, see equation 35, page
    #: 1021
    DEFINED_FOR_STANDARD_DEVIATION_TYPES = set([
        const.StdDev.TOTAL,
        const.StdDev.INTER_EVENT,
        const.StdDev.INTRA_EVENT
    ])

    #: No site parameters are needed
    REQUIRES_SITES_PARAMETERS = set()

    #: Required rupture parameter is only magnitude, see equation 30 page
    #: 1021.
    REQUIRES_RUPTURE_PARAMETERS = set(('mag', ))

    #: Required distance measure is closest distance to rupture, see equation
    #: 30 page 1021.
    REQUIRES_DISTANCES = set(('rjb', ))

    def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        """
        See :meth:`superclass method
        <.base.GroundShakingIntensityModel.get_mean_and_stddevs>`
        for spec of input and result values.
        """
        #assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
        #           for stddev_type in stddev_types)

        C = self.COEFFS[imt]
        mean = self._compute_mean(C, rup, dists.rjb)
        if isinstance(imt, (SA, PGA)): # Convert from m/s**2 to g
            mean -= np.log(g)
        elif isinstance(imt, PGV): # Convert from m/s to cm/s
            mean += np.log(100.0)
        stddevs = self._get_stddevs(C, stddev_types, rup.mag,
                                    dists.rjb.shape)

        return mean, stddevs

    def _compute_mean(self, C, rup, rjb):
        """
        Compute mean value according to equation 30, page 1021.
        """
        mean = (C['c1'] +
                self._compute_magnitude_term(C, rup) +
                self._compute_distance_term(C, rup, rjb))
        return mean

    def _get_stddevs(self, C, stddev_types, mag, num_sites):
        """
        Return total standard deviation as for equation 35, page 1021.
        """
        stddevs = []
        for stddev_type in stddev_types:
            assert stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
            if stddev_type == const.StdDev.TOTAL:
                sigma_t = np.sqrt(C['sigma'] ** 2. + C['tau'] ** 2.)
                stddevs.append(sigma_t + np.zeros(num_sites))
            elif stddev_type == const.StdDev.INTRA_EVENT:
                stddevs.append(C['sigma'] + np.zeros(num_sites))
            elif stddev_type == const.StdDev.INTER_EVENT:
                stddevs.append(C['tau'] + np.zeros(num_sites))
        return stddevs

    def _compute_magnitude_term(self, C, rup):
        """
        This computes the term f1 equation 8 Drouet & Cotton (2015)
        """
        return C['c2'] * (rup.mag - 8.0) + C['c3'] * (rup.mag - 8.0) ** 2

    def _compute_distance_term(self, C, rup, rjb):
        """
        This computes the term f2 equation 8 Drouet & Cotton (2015)
        """
        return (C['c4'] + C['c5'] * rup.mag) * np.log(
            np.sqrt(rjb ** 2. + C['c6'] ** 2.)) + C['c7'] * rjb

    #: Coefficient tables are constructed from the electronic suplements of
    #: the original paper.
    COEFFS = CoeffsTable(sa_damping=5, table="""\
    IMT          c1         c2         c3         c4        c5        c6         c7     sigma       tau
    pgv    0.457169   0.230048  -0.140795  -1.785544  0.162408  7.201555  -0.002832  0.526166  0.393025
    pga    2.575109  -0.243140  -0.155164  -1.807995  0.156084  7.629410  -0.003996  0.625075  0.481226
    0.010  2.586201  -0.242271  -0.154703  -1.808695  0.156141  7.623170  -0.004020  0.630187  0.481707
    0.020  2.717172  -0.261739  -0.152580  -1.859646  0.160136  7.773640  -0.004048  0.688591  0.486003
    0.030  2.930319  -0.236637  -0.147729  -1.831982  0.155727  7.918015  -0.004558  0.751663  0.490319
    0.040  3.156209  -0.266436  -0.143032  -1.915769  0.162181  8.224381  -0.004350  0.779379  0.498289
    0.050  3.268731  -0.189794  -0.140262  -1.753521  0.144779  7.767181  -0.004910  0.780218  0.497182
    0.075  3.352582  -0.164010  -0.140254  -1.691894  0.140095  7.428414  -0.004756  0.706350  0.493511
    0.100  3.455122  -0.203575  -0.148680  -1.708867  0.139700  7.707583  -0.004261  0.652456  0.490010
    0.150  3.456514  -0.169395  -0.160434  -1.607720  0.131021  7.274064  -0.004025  0.587130  0.480912
    0.200  3.480893  -0.155262  -0.168476  -1.646459  0.132556  7.424609  -0.002871  0.556933  0.462619
    0.250  3.358985  -0.255601  -0.194574  -1.669187  0.134462  7.753731  -0.002732  0.533650  0.458696
    0.300  3.115954  -0.237559  -0.215762  -1.451276  0.114182  7.212529  -0.003761  0.526336  0.452876
    0.400  2.806835  -0.340296  -0.250121  -1.418762  0.114054  6.837724  -0.004081  0.519411  0.440790
    0.500  2.837393  -0.355473  -0.271003  -1.453916  0.111753  7.298391  -0.003037  0.512892  0.427910
    0.750  2.383076  -0.374649  -0.298428  -1.472297  0.117984  7.051676  -0.002899  0.507442  0.405868
    1.000  2.070536  -0.263869  -0.303220  -1.410898  0.117144  6.815268  -0.003307  0.511352  0.384417
    1.250  1.944386  -0.196142  -0.309115  -1.408815  0.116519  6.904435  -0.003017  0.511909  0.376152
    1.500  1.973072  -0.160616  -0.313180  -1.493457  0.122469  7.427893  -0.002316  0.511871  0.370833
    1.750  1.747518  -0.129961  -0.320672  -1.400692  0.116855  7.143261  -0.003402  0.508641  0.361738
    2.000  1.667278  -0.083863  -0.319818  -1.405853  0.114769  7.128404  -0.003174  0.505025  0.353357
    3.000  1.292331   0.312316  -0.263539  -1.464213  0.130085  6.416692  -0.002621  0.512370  0.344082
    """)

class DrouetBrazil2015withDepth(DrouetBrazil2015):
    """
    Implements GMPE developed by S. Drouet unpublished for Brazil based on the
    method described in Douet & Cotton (2015) BSSA doi: 10.1785/0120140240.
    Model with magnitude-dependent depth distribution and depth-dependent
    stress distribution
    """
    #: Required rupture parameter is only magnitude, see equation 30 page
    #: 1021.
    REQUIRES_RUPTURE_PARAMETERS = set(('mag', 'hypo_depth'))

    def _compute_distance_term(self, C, rup, rjb):
        """
        This computes the term f2 equation 8 Drouet & Cotton (2015)
        """
        return (C['c4'] + C['c5'] * rup.mag) * np.log(
            np.sqrt(rjb ** 2. + rup.hypo_depth ** 2.)) + C['c6'] * rjb

    #: Coefficient tables are constructed from the electronic supplements of
    #: the original paper.
    COEFFS = CoeffsTable(sa_damping=5, table="""\
    IMT            c1        c2        c3        c4       c5        c6    sigma      tau
    pgv      1.296890  0.307642 -0.155927 -1.727495 0.128331 -0.002079 0.507600 0.474106
    pga      3.393123 -0.139286 -0.169796 -1.695729 0.116686 -0.003491 0.609317 0.573384
    0.010000 3.405471 -0.138584 -0.169352 -1.696947 0.116766 -0.003513 0.614578 0.574138
    0.020000 3.512058 -0.152120 -0.167461 -1.725659 0.119114 -0.003688 0.674715 0.580870
    0.030000 3.695045 -0.123957 -0.162204 -1.688566 0.114624 -0.004259 0.738161 0.585416
    0.040000 3.856703 -0.140850 -0.157819 -1.721169 0.117496 -0.004400 0.768733 0.597967
    0.050000 4.085005 -0.080044 -0.155124 -1.628430 0.103907 -0.004434 0.766573 0.593676
    0.075000 4.242841 -0.056117 -0.154044 -1.601802 0.100823 -0.004025 0.691819 0.587975
    0.100000 4.278594 -0.087360 -0.162257 -1.582776 0.098667 -0.003803 0.638952 0.582946
    0.150000 4.358265 -0.063108 -0.173452 -1.530628 0.092988 -0.003244 0.572059 0.565676
    0.200000 4.385790 -0.041305 -0.181717 -1.553364 0.092537 -0.002134 0.539692 0.543338
    0.250000 4.162408 -0.135372 -0.207677 -1.531599 0.092853 -0.002369 0.520273 0.534858
    0.300000 3.991746 -0.134064 -0.227479 -1.386035 0.078380 -0.002944 0.511192 0.519484
    0.400000 3.746791 -0.248433 -0.261965 -1.396553 0.080903 -0.002965 0.502079 0.499394
    0.500000 3.728968 -0.246365 -0.282618 -1.384885 0.074963 -0.002185 0.495203 0.483281
    0.750000 3.298203 -0.278081 -0.310518 -1.431361 0.083594 -0.001906 0.488881 0.452970
    1.000000 2.966504 -0.185734 -0.315139 -1.405592 0.087333 -0.002169 0.492014 0.426399
    1.250000 2.810007 -0.117505 -0.320881 -1.394705 0.086834 -0.001981 0.493565 0.416078
    1.500000 2.760153 -0.068911 -0.325231 -1.426076 0.089344 -0.001624 0.494063 0.410019
    1.750000 2.538491 -0.052014 -0.332247 -1.367936 0.087458 -0.002541 0.491032 0.398572
    2.000000 2.485648 -0.001655 -0.331560 -1.370772 0.084286 -0.002297 0.487175 0.392599
    3.000000 2.205128  0.365128 -0.276508 -1.504904 0.104957 -0.001304 0.492143 0.384311
    """)

max_stars_repo_name: gfzriesgos/shakyground-lfs

max_issues_repo_name: gfzriesgos/shakyground-lfs

max_forks_repo_name: gfzriesgos/shakyground-lfs

reason: import numpy,from scipy


text:
def readline():
    import sys

    return sys.stdin.buffer.readline().rstrip()

def readline_ints():
    (*ints,) = map(
        int,
        readline().split(),
    )
    return ints

import numpy as np

def get_result(n, a):
    b = list(range(n))
    for i in a[::-1]:
        b[i - 1], b[i] = b[i], b[i - 1]
    return b

def solve(n, a, d):
    a = get_result(n, a)
    a = np.array(a)
    res = f(a, d) + 1
    print(*res, sep="\n")

def f(a, n):
    if n == 0:
        return np.arange(a.size)
    b = f(a, n // 2)
    b = b[b]
    if n & 1:
        b = a[b]
    return b

def main():
    n, _, d = readline_ints()
    a = readline_ints()
    solve(n, a, d)

if __name__ == "__main__":
    main()

max_stars_repo_name: kagemeka/atcoder-submissions

max_issues_repo_name: kagemeka/atcoder-submissions

max_forks_repo_name: kagemeka/atcoder-submissions


text:
#Program simulate queueing system in parallel
#Copyright (C) 2016 by Araik Tamazian, Viet Duc Nguyen

import pandas as pd
import numpy as np
from gg1_function import simulate_gg1, qexp_rate, rand_qexp, simulate_mm1
import matplotlib.pyplot as plt
import timeit
import multiprocessing
from joblib import Parallel, delayed

#Import external file containing empirical data (time and size)
infile = 'ses_20081013.txt'
data = pd.read_csv(infile,delim_whitespace = True, header=None, na_filter = True)      #Read file without header and using space ' ' to seperate between columns
data.columns = ['time','ip', 'numcon', 'size']   # Sign names for columns
data = data[data.time >= 0]                    #Drop negative time variable
time1 = np.asarray(data['time'])
time1 = np.diff(time1)
time = np.insert(time1,0,0.0)
time  = time/np.mean(time)
ssize = data['size']
ssize = ssize/np.mean(ssize)
print 'len time:', len(time)
print 'len size: ', len(ssize)

#Import surrogate time
timess = pd.read_csv('sur_time_20081013_2.txt')
timess = np.asarray(timess)
timess = timess.flatten()
ssizess = pd.read_csv('sur_size_20081013_2.txt')
ssizess = np.asarray(ssizess)
ssizess = ssizess.flatten()
print 'len timess:',len(timess)
print 'len ssizess: ',len(ssizess)

#Declare variables
res_mm1 =[]                                                     #Variale gets results of simulating MM1
n = 10**6                                                       #Number of data samples (data lines) needed simulating
c=np.logspace(np.log10(1.2), np.log10(10),10)                   #Values of throughput
inputs = range(len(c))                                          #Variable for loop in parallel

#Function for simulating queueing system with empirical data in parallel via throughput
def processInput(i):
    return simulate_gg1(n,time,ssize/c[i])

#Function for simulating queueing system with surrogate data in parallel via throughput
def processInputs(i):
    return simulate_gg1(n,timess,ssizess/c[i])
if __name__ == '__main__':
    start_time = timeit.default_timer()                                         #Variale for calculate time of simulation

    #Simulate queueing system with empirical data, surrogate data and MM1
    res_real  = Parallel(n_jobs=3)(delayed(processInput)(i) for i in inputs)    #With empirical data
    print 'empirical data:', res_real

    res_sur  = Parallel(n_jobs=3)(delayed(processInputs)(i) for i in inputs)    #With surrogate data
    print 'Surrogate data:', res_sur

    for i in range(len(c)):                                                     #Queueing system MM1
        res_mm1.append(simulate_mm1(1/np.mean(time), 1/np.mean(ssize/c[i])))

    #Convert results from array into data.frame using thoughput like index and save into txt files
    res_real = pd.DataFrame(data = res_real,index=c, columns=['Ur','Wr','Lr'])
    res_sur = pd.DataFrame(data = res_sur,index=c, columns=['Us','Ws','Ls'])
    res_mm1 = pd.DataFrame(data = res_mm1,index=c, columns=['Um','Wm','Lm'])

    res_real.to_csv("simres_real.txt", sep=" ", header=False)
    res_sur.to_csv("simres_sim.txt", sep=" ", header=False)
    res_mm1.to_csv("simres_mm1.txt", sep=" ", header=False)

    #------Plot multigraphs in one panel (subplot)-------
    #Plot Utilization versus throughput
    plt.subplot(1,3,1)                                              #plt.subplot(row,column, position)
    plt.plot(c,res_real['Ur'], color = 'black',linewidth = 2)
    plt.plot(c,res_sur['Us'], color = 'red',linewidth = 2)
    plt.plot(c,res_mm1['Um'], color = 'green',linewidth = 2)
    plt.xlabel('Channel throughput, C', size = 15)
    plt.ylabel('Utilization, U', size = 15)
    plt.yscale('log')
    plt.xscale('log')
    plt.xticks(color='black', size=15)
    plt.yticks(color='black', size=15)
    plt.legend(loc = "upper right")

    #Plot Average sojourn time in system versus throughput
    plt.subplot(1,3,2)
    plt.plot(c,res_real['Wr'], color = 'black',linewidth = 2)
    plt.plot(c,res_sur['Ws'], color = 'red',linewidth = 2)
    plt.plot(c,res_mm1['Wm'], color = 'green',linewidth = 2)
    plt.xlabel('Channel throughput, C', size = 15)
    plt.ylabel('Average soujourn time,W', size = 15)
    plt.yscale('log')
    plt.xscale('log')
    plt.xticks(color='black', size=15)
    plt.yticks(color='black', size=15)
    plt.legend(loc = "upper right")

    #Plot Average number of requests in system versus throughput
    plt.subplot(1,3,3)
    plt.loglog(c,res_real['Lr'], color = 'black',linewidth = 2)
    plt.loglog(c,res_sur['Ls'], color = 'red',linewidth = 2)
    plt.loglog(c,res_mm1['Lm'], color = 'green',linewidth = 2)
    plt.xlabel('Channel throughput, C', size = 15)
    plt.ylabel('Average request in system, L', size = 15)
    plt.xticks(color='black', size=15)
    plt.yticks(color='black', size=15)
    plt.yscale('log')
    plt.xscale('log')
    plt.legend(loc = "upper right")

    stop_time = timeit.default_timer()
    time_simulation = stop_time - start_time
    print 'simulation time: ', time_simulation

    plt.savefig('n206.png')
    plt.show()

#---------------end-----------------------------------------------------------------------------------------------------


text:
import numpy as np

from ...misc import *
from ....types import *
from ..base import Converter

@Converter.register("string")
class StringConverter(Converter):
    def _fit(self) -> "Converter":
        assert self.info.is_valid and self.info.is_string and self.info.need_transform
        self._transform_dict = self._recognizer.transform_dict
        reversed_dict = {v: k for k, v in self._transform_dict.items()}
        self._reverse_transform_dict = reversed_dict
        return self

    def _convert(self, flat_arr: flat_arr_type) -> np.ndarray:
        converted = [self._transform_dict.get(elem, 0) for elem in flat_arr]
        return np.array(converted, np_float_type)

    def _recover(self, flat_arr: flat_arr_type) -> np.ndarray:
        recovered = [self._reverse_transform_dict[elem] for elem in flat_arr]
        return np.array(recovered, np.str)

__all__ = ["StringConverter"]