text:
import pathlib
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from scipy import stats

ops = ['a', 'd', 're', 'r', 's', 'sum', 'su', 'pr']

#Currently we read all the files for each plot, this is kinda slow and could be fixed.

def optoprettytitle(op):
    if op == 'a':
        return 'Add value'
    if op == 'd':
        return 'Delete index'
    if op == 're':
        return 'Remove value'
    if op == 'r':
        return 'Rank value'
    if op == 's':
        return 'Select index'
    if op == 'sum':
        return 'Sum values'
    if op == 'su':
        return 'Successor of value'
    if op == 'pr':
        return 'Predecessor of value'

def createplotforoperation(op):
    def mad(data, axis=None):
        return np.mean(np.absolute(data - np.mean(data, axis)), axis)

    def dfcsv(name):
        data = pd.read_csv(f"../results/timings/{name}.large.csv")
        df = pd.DataFrame(data[data["Operation"] == op],columns=['Time(s)','Size'])
        df['median'] = df['Time(s)'].rolling(1000).median()
        # Below line removes outliers but is not necessary
        #df['std'] = df['Time(s)'].rolling(10).std()
        #df = df[(df['Time(s)'] <= df['median']+0.5*df['std']) & (df['Time(s)'] >= df['median']-3*df['std'])]
        return df

    def addplot(name, ax):
        df = dfcsv(name)
        df.plot(x ='Size', y='median', kind = 'line', ax = ax, label=name)

    df = dfcsv("SortedContainer")
    ax = df.plot(x ='Size', y='median', kind = 'line', label="SortedContainer")

    addplot("SortedArray", ax)
    addplot("Blist", ax)
    addplot("AutoLoad", ax)
    #addplot("Skiplist", ax)
    #addplot("RBSTree", ax)
    #addplot("SortedCollection", ax)

    locs,labels = plt.xticks()
    plt.xticks(locs, map(lambda x: "%i" % x, locs))
    plt.yticks(locs, map(lambda x: "%f" % x, locs))
    plt.xlabel("N")
    plt.ylabel("Time(ns)")
    plt.grid(True, which="both", linestyle='--')
    ax.set_xscale('log')
    ax.set_yscale('log')

    # Shrink current axis by 20%
    box = ax.get_position()
    ax.set_position([box.x0, box.y0, box.width * 0.8, box.height])

    # Put a legend to the right of the current axis
    ax.legend(loc='center left', bbox_to_anchor=(1, 0.5))

    # Resize plot figure
    plt.gcf().set_size_inches(10, 6)

    plt.title(f'Median runtimes for operation "{optoprettytitle(op)}" over N')
    pathlib.Path('../results/graphs/').mkdir(parents=True, exist_ok=True)
    plt.savefig(f'../results/graphs/{op}.png')

if __name__ == '__main__':
    for op in ops:
        createplotforoperation(op)

max_stars_repo_name: Jakobis/OrderedSequences

max_issues_repo_name: Jakobis/OrderedSequences

max_forks_repo_name: Jakobis/OrderedSequences

reason: import numpy,from scipy


text:
import torch
import torch.nn as nn
import torch.nn.functional as F
import numpy as np

device = torch.device("cuda" if torch.cuda.is_available() else "cpu")

class Agent(nn.Module):
    def __init__(self, **kwargs):
        super(Agent, self).__init__(**kwargs)

        self.init_args = kwargs

        self.h = torch.zeros(1, 256)

        self.norm1 = nn.BatchNorm2d(4)
        self.norm2 = nn.BatchNorm2d(32)

        self.conv1 = nn.Conv2d(4, 32, 4, stride=2, padding=1)
        self.conv2 = nn.Conv2d(32, 32, 3, stride=2, padding=1)
        self.conv3 = nn.Conv2d(32, 32, 3, stride=2, padding=1)
        self.conv4 = nn.Conv2d(32, 32, 3, stride=2, padding=1)

        self.gru = nn.GRUCell(1152, 256)
        self.policy = nn.Linear(256, 4)
        self.value = nn.Linear(256, 1)

        self.prelu1 = nn.PReLU()
        self.prelu2 = nn.PReLU()
        self.prelu3 = nn.PReLU()
        self.prelu4 = nn.PReLU()

        nn.init.xavier_uniform_(self.conv1.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.conv1.bias, 0.01)
        nn.init.xavier_uniform_(self.conv2.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.conv2.bias, 0.01)
        nn.init.xavier_uniform_(self.conv3.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.conv3.bias, 0.01)
        nn.init.xavier_uniform_(self.conv4.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.conv4.bias, 0.01)
        nn.init.constant_(self.gru.bias_ih, 0)
        nn.init.constant_(self.gru.bias_hh, 0)
        nn.init.xavier_uniform_(self.policy.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.policy.bias, 0.01)
        nn.init.xavier_uniform_(self.value.weight, gain=nn.init.calculate_gain('leaky_relu'))
        nn.init.constant_(self.value.bias, 0.01)

        self.train()

    def reset(self):
        self.h = torch.zeros(1, 256)

    def clone(self, num=1):
        return [ self.clone_one() for _ in range(num) ]

    def clone_one(self):
        return Agent(**self.init_args)

    def forward(self, state):
        state = state.view(1, 4, 96, 96)
        state = self.norm1(state)

        data = self.prelu1(self.conv1(state))
        data = self.prelu2(self.conv2(data))
        data = self.prelu3(self.conv3(data))
        data = self.prelu4(self.conv4(data))

        data = self.norm2(data)
        data = data.view(1, -1)

        h = self.gru(data, self.h)
        self.h = h.detach()

        pre_policy = h.view(-1)

        policy = F.softmax(self.policy(pre_policy))
        value = self.value(pre_policy)

        return policy, value

ext: py

reason: import numpy


text:
import sys

sys.path.append('')
import pytorch_lightning as pl
import torch
import torch.nn as nn
from typing import *
import numpy as np
import random

from languasito.utils import Encodings, mask_concat
from languasito.modules import WordGram, LinearNorm, CosineLoss, WordDecoder

class Languasito(pl.LightningModule):
    def __init__(self, encodings: Encodings):
        super().__init__()
        NUM_FILTERS = 512
        RNN_SIZE = 256
        CHAR_EMB_SIZE = 128
        ATT_DIM = 64
        NUM_HEADS = 8

        self._wg = WordGram(len(encodings.char2int), num_langs=1, num_filters=512, num_layers=5)
        self._rnn_fw = nn.LSTM(NUM_FILTERS // 2, RNN_SIZE, num_layers=3, batch_first=True, bidirectional=False)
        self._rnn_bw = nn.LSTM(NUM_FILTERS // 2, RNN_SIZE, num_layers=3, batch_first=True, bidirectional=False)
        self._linear_out = LinearNorm(RNN_SIZE * 2, NUM_FILTERS // 2)
        self._early_stop_meta_val = 0
        self._res = {"b_loss": 9999}
        self._start_stop = nn.Embedding(2, NUM_FILTERS // 2)
        self._epoch_results = None
        self._loss_function = nn.CrossEntropyLoss(ignore_index=0)
        self._repr1_ff = nn.Sequential(nn.LayerNorm(RNN_SIZE), nn.Linear(RNN_SIZE, NUM_FILTERS), nn.ReLU(),
                                       nn.LayerNorm(NUM_FILTERS), nn.Linear(NUM_FILTERS, NUM_FILTERS), nn.ReLU())
        self._repr2_ff = nn.Sequential(nn.LayerNorm(ATT_DIM * NUM_HEADS), nn.Linear(ATT_DIM * NUM_HEADS, NUM_FILTERS),
                                       nn.ReLU(),
                                       nn.LayerNorm(NUM_FILTERS), nn.Linear(NUM_FILTERS, NUM_FILTERS), nn.ReLU())
        self._key = nn.Sequential(nn.Linear(RNN_SIZE, ATT_DIM), nn.Tanh())
        self._value = nn.Sequential(nn.Linear(RNN_SIZE, ATT_DIM), nn.Tanh())
        self._att_fn_fw = nn.MultiheadAttention(RNN_SIZE, NUM_HEADS, kdim=ATT_DIM, vdim=ATT_DIM)
        self._att_fn_bw = nn.MultiheadAttention(RNN_SIZE, NUM_HEADS, kdim=ATT_DIM, vdim=ATT_DIM)
        cond_size = NUM_FILTERS
        self._word_reconstruct = WordDecoder(cond_size, CHAR_EMB_SIZE, len(encodings.word_decomposer._tok2int) + 4)
        self._cosine_loss = CosineLoss()

    def forward(self, X, return_w=False, imagine=False):
        x_words_chars = X['x_word_char']
        x_words_case = X['x_word_case']
        x_lang_word = X['x_lang_word']
        x_sent_len = X['x_sent_len']
        x_word_len = X['x_word_len']
        x_word_masks = X['x_word_masks']
        x_max_len = X['x_max_len']
        char_emb_packed = self._wg(x_words_chars, x_words_case, x_lang_word, x_word_masks, x_word_len)

        blist_char = []

        sl = x_sent_len.cpu().numpy()
        pos = 0
        for ii in range(x_sent_len.shape[0]):
            slist_char = []
            slist_char.append(
                self._start_stop(torch.zeros((1), dtype=torch.long, device=self._get_device())))
            for jj in range(sl[ii]):
                slist_char.append(char_emb_packed[pos, :].unsqueeze(0))
                pos += 1

            slist_char.append(
                self._start_stop(torch.ones((1), dtype=torch.long, device=self._get_device())))

            for jj in range(x_max_len - sl[ii]):
                slist_char.append(torch.zeros((1, 512 // 2),
                                              device=self._get_device(), dtype=torch.float))

            sent_emb = torch.cat(slist_char, dim=0)
            blist_char.append(sent_emb.unsqueeze(0))

        char_emb = torch.cat(blist_char, dim=0)
        out_fw, _ = self._rnn_fw(char_emb)
        out_bw, _ = self._rnn_bw(torch.flip(char_emb, [1]))
        out_bw = torch.flip(out_bw, [1])
        lexical = char_emb[:, 1:-1, :]
        out_fw = out_fw + char_emb
        out_bw = out_bw + char_emb
        pre_context = torch.cat([out_fw[:, :-2, :], out_bw[:, 2:, :]], dim=-1)
        context = torch.tanh(self._linear_out(pre_context))

        concat = torch.cat([lexical, context, pre_context], dim=-1)

        y = {'lexical': lexical, 'context': context, 'emb': concat}  # , 'sent': sent}

        if return_w:
            att_value = self._apply_masked_attention(out_fw[:, :-2, :], out_bw[:, 2:, :])
            # att_value = torch.zeros_like(att_value)
            # context = torch.zeros_like(context)
            # repr1 = self._repr1_ff(context)
            repr2 = self._repr2_ff(att_value)
            # cond = torch.cat([repr1, repr2], dim=-1)
            # cond = mask_concat([repr1, repr2], 0.33, self.training, self._get_device())
            cond = repr2

            # cond = repr1
            cond_packed = []
            for ii in range(x_sent_len.shape[0]):
                for jj in range(x_sent_len[ii]):
                    cond_packed.append(cond[ii, jj].unsqueeze(0))
            cond_packed = torch.cat(cond_packed, dim=0)
            if imagine:
                x_char_pred = self._word_reconstruct(cond_packed, gs_chars=None)
            else:
                x_char_pred = self._word_reconstruct(cond_packed, gs_chars=X['x_word_targets'])
            y['x_char_pred'] = x_char_pred

        return y

    def _apply_masked_attention(self, fw, bw):
        # forward
        att_query = fw
        att_mask = np.ones((fw.shape[1], fw.shape[1]), dtype=np.float)
        for ii in range(fw.shape[1]):
            for jj in range(ii + 1, fw.shape[1]):
                att_mask[ii, jj] = 0
        att_mask = torch.tensor(att_mask, device=self._get_device())
        att_key = self._key(fw)
        att_val = self._value(fw)
        att_mask = att_mask.float().masked_fill(att_mask == 0, float('-inf')).masked_fill(att_mask == 1, float(0.0))
        att_value_fw, _ = self._att_fn_fw(att_query.permute(1, 0, 2), att_key.permute(1, 0, 2),
                                          att_val.permute(1, 0, 2), attn_mask=att_mask)

        att_query = bw
        att_mask = np.ones((bw.shape[1], bw.shape[1]), dtype=np.float)
        for ii in range(bw.shape[1]):
            for jj in range(0, ii):
                att_mask[ii, jj] = 0
        att_mask = torch.tensor(att_mask, device=self._get_device())
        att_mask = att_mask.float().masked_fill(att_mask == 0, float('-inf')).masked_fill(att_mask == 1, float(0.0))
        att_key = self._key(bw)
        att_val = self._value(bw)

        att_value_bw, _ = self._att_fn_bw(att_query.permute(1, 0, 2), att_key.permute(1, 0, 2),
                                          att_val.permute(1, 0, 2), attn_mask=att_mask)

        return torch.cat([att_value_fw, att_value_bw], dim=-1).permute(1, 0, 2)

    def training_step(self, batch, batch_idx):
        Y = self.forward(batch, return_w=True)
        x_char_target = batch['x_word_targets'][:, 1:]
        x_char_pred = Y['x_char_pred']
        loss_rec = self._loss_function(x_char_pred.reshape(-1, x_char_pred.shape[2]), x_char_target.reshape(-1))
        # y_lexical = Y['lexical']
        # y_context = Y['context']
        # sl = batch['x_sent_len'].detach().cpu().numpy()
        # word_repr = []
        # sent_repr = y_sent
        # for ii in range(sl.shape[0]):
        #     for jj in range(sl[ii]):
        #         if True:  # random.random() < 0.15:
        #             word_repr.append(y_lexical[ii, jj].unsqueeze(0))
        #             word_repr.append(y_context[ii, jj].unsqueeze(0))
        #
        # word_repr = torch.cat(word_repr, dim=0)
        # word_repr = word_repr.reshape(-1, 2, word_repr.shape[1])
        # loss_cosine = self._cosine_loss(word_repr)
        #
        # # sent_repr = sent_repr.reshape(-1, 2, sent_repr.shape[1])
        # # loss_sent = self._ge2e_sent(sent_repr)
        return loss_rec  # + loss_cosine

    def validation_step(self, batch, batch_idx):
        Y = self.forward(batch, return_w=True)
        x_char_target = batch['x_word_targets'][:, 1:]
        x_char_pred = Y['x_char_pred']

        loss_rec = self._loss_function(x_char_pred.reshape(-1, x_char_pred.shape[2]), x_char_target.reshape(-1))

        # y_lexical = Y['lexical']
        # y_context = Y['context']
        # sl = batch['x_sent_len'].detach().cpu().numpy()
        # word_repr = []
        # # sent_repr = y_sent
        # for ii in range(sl.shape[0]):
        #     for jj in range(sl[ii]):
        #         if True:  # random.random() < 0.15:
        #             word_repr.append(y_lexical[ii, jj].unsqueeze(0))
        #             word_repr.append(y_context[ii, jj].unsqueeze(0))
        #
        # word_repr = torch.cat(word_repr, dim=0)
        # word_repr = word_repr.reshape(-1, 2, word_repr.shape[1])
        # loss_cosine = self._cosine_loss(word_repr)
        # return {'total_loss': loss_rec + loss_cosine}
        return {'total_loss': loss_rec}

    def validation_epoch_end(self, outputs: List[Any]) -> None:

        loss = 0
        for output in outputs:
            loss += output['total_loss']
        loss /= len(outputs)

        res = {'val_loss': loss}
        self._epoch_results = self._compute_early_stop(res)
        self.log('val/early_meta', self._early_stop_meta_val)
        self.log('val/loss', loss)

    def configure_optimizers(self):
        return torch.optim.AdamW(self.parameters())

    def _compute_early_stop(self, res):
        if res["val_loss"] < self._res['b_loss']:
            self._early_stop_meta_val += 1
            self._res['b_loss'] = res["val_loss"]
            res['best_loss'] = True
        return res

    def _get_device(self):
        if self._start_stop.weight.device.type == 'cpu':
            return 'cpu'
        return '{0}:{1}'.format(self._start_stop.weight.device.type, str(self._start_stop.weight.device.index))

    def load(self, filename: str):
        self.load_state_dict(torch.load(filename, map_location='cpu'))

max_stars_repo_name: AliOsamaHassan/NLP-Cube

max_issues_repo_name: davidegazze/NLP-Cube

max_forks_repo_name: davidegazze/NLP-Cube


text:
import math
import pkg_resources
import random

# 3rd party modules
from gym import spaces
import gym
import numpy as np

def get_chance(x):
    e = math.exp(1)
    return (1.0 + e) / (1. + math.exp(x + 1))

class BananaEnv(gym.Env):
    def __init__(self):
        # General variables defining the environment
        self.MAX_PRICE = 2.0
        self.TOTAL_TIME_STEPS = 2

        self.curr_step = -1
        self.is_banana_sold = False

        # Define what the agent can do
        # Sell at 0.00 EUR, 0.10 Euro, ..., 2.00 Euro
        self.action_space = spaces.Discrete(21)

        # Observation is the remaining time
        low = np.array([
            0.0,  # remaining_tries
        ])
        high = np.array([
            self.TOTAL_TIME_STEPS,  # remaining_tries
        ])
        self.observation_space = spaces.Box(low, high, dtype=np.float32)

        # Store what the agent tried
        self.curr_episode = -1
        self.action_episode_memory = []

    def step(self, action):
        """
        The agent takes a step in the environment.
        Parameters
        ----------
        action : int
        Returns
        -------
        ob, reward, episode_over, info : tuple
            ob (object) :
                an environment-specific object representing your observation of
                the environment.
            reward (float) :
                amount of reward achieved by the previous action. The scale
                varies between environments, but the goal is always to increase
                your total reward.
            episode_over (bool) :
                whether it's time to reset the environment again. Most (but not
                all) tasks are divided up into well-defined episodes, and done
                being True indicates the episode has terminated. (For example,
                perhaps the pole tipped too far, or you lost your last life.)
            info (dict) :
                 diagnostic information useful for debugging. It can sometimes
                 be useful for learning (for example, it might contain the raw
                 probabilities behind the environment's last state change).
                 However, official evaluations of your agent are not allowed to
                 use this for learning.
        """
        if self.is_banana_sold:
            raise RuntimeError("Episode is done")
        self.curr_step += 1
        self._take_action(action)
        reward = self._get_reward()
        ob = self._get_state()
        return ob, reward, self.is_banana_sold, {}

    def _take_action(self, action):
        self.action_episode_memory[self.curr_episode].append(action)
        self.price = ((float(self.MAX_PRICE) / (self.action_space.n - 1)) *
                      action)

        chance_to_take = get_chance(self.price)
        banana_is_sold = (random.random() < chance_to_take)

        if banana_is_sold:
            self.is_banana_sold = True

        remaining_steps = self.TOTAL_TIME_STEPS - self.curr_step
        time_is_over = (remaining_steps <= 0)
        throw_away = time_is_over and not self.is_banana_sold
        if throw_away:
            self.is_banana_sold = True  # abuse this a bit
            self.price = 0.0

    def _get_reward(self):
        """Reward is given for a sold banana."""
        if self.is_banana_sold:
            return self.price - 1
        else:
            return 0.0

    def reset(self):
        """
        Reset the state of the environment and returns an initial observation.
        Returns
        -------
        observation (object): the initial observation of the space.
        """
        self.curr_step = -1
        self.curr_episode += 1
        self.action_episode_memory.append([])
        self.is_banana_sold = False
        self.price = 1.00
        return self._get_state()

    def _render(self, mode='human', close=False):
        return

    def _get_state(self):
        """Get the observation."""
        ob = [self.TOTAL_TIME_STEPS - self.curr_step]
        return ob

    def seed(self, seed):
        random.seed(seed)
        np.random.seed


text:
import tensorflow as tf
import numpy as np
from PIL import Image
import pickle
import matplotlib.pyplot as plt
import cv2

import math

def parse_output(heatmap_data,offset_data, threshold):
  '''
  Input:
    heatmap_data - hetmaps for an image. Three dimension array
    offset_data - offset vectors for an image. Three dimension array
    threshold - probability threshold for the keypoints. Scalar value
  Output:
    array with coordinates of the keypoints and flags for those that have
    low probability
  '''
  joint_num = heatmap_data.shape[-1]
  pose_kps = np.zeros((joint_num,3), np.uint32)

  for i in range(heatmap_data.shape[-1]):
      joint_heatmap = heatmap_data[...,i]
      max_val_pos = np.squeeze(np.argwhere(joint_heatmap==np.max(joint_heatmap)))
      remap_pos = np.array(max_val_pos/8*257,dtype=np.int32)
      pose_kps[i,0] = int(remap_pos[0] + offset_data[max_val_pos[0],max_val_pos[1],i])
      pose_kps[i,1] = int(remap_pos[1] + offset_data[max_val_pos[0],max_val_pos[1],i+joint_num])
      max_prob = np.max(joint_heatmap)

      if max_prob > threshold:
        if pose_kps[i,0] < 257 and pose_kps[i,1] < 257:
          pose_kps[i,2] = 1
  return pose_kps

def new_coordinates_after_resize_img(original_size, new_size, original_coordinate):
  original_size = np.array(original_size)
  new_size = np.array(new_size)
  original_coordinate = np.array(original_coordinate)
  xy = original_coordinate/(original_size/new_size)
  x, y = int(xy[0]), int(xy[1])
  return (x, y)

def poseNet(input):
    model_path = "pretrained_model/posenet_mobilenet_v1_100_257x257_multi_kpt_stripped.tflite"
    interpreter = tf.lite.Interpreter(model_path=model_path)
    interpreter.allocate_tensors()
    input_details = interpreter.get_input_details()
    output_details = interpreter.get_output_details()
    input_shape = input_details[0]['shape']

    height,width = input.shape[0],input.shape[1]
    test_img = cv2.resize(input,(input_shape[1],input_shape[2]))
    h,w = test_img.shape[0],test_img.shape[1]
    template_input = np.expand_dims(test_img.copy(), axis=0)
    floating_model = input_details[0]['dtype'] == np.float32

    if floating_model:
        template_input = (np.float32(template_input) - 127.5) / 127.5

    interpreter.set_tensor(input_details[0]['index'], template_input)
    interpreter.invoke()
    template_output_data = interpreter.get_tensor(output_details[0]['index'])
    template_offset_data = interpreter.get_tensor(output_details[1]['index'])
    template_heatmaps = np.squeeze(template_output_data)
    template_offsets = np.squeeze(template_offset_data)
    template_show = np.squeeze((template_input.copy()*127.5+127.5)/255.0)
    template_show = np.array(template_show*255,np.uint8)
    template_kps = parse_output(template_heatmaps,template_offsets,0.3)

    left_neck_x = template_kps[6][1]
    left_neck_y = template_kps[6][0]

    right_neck_x = template_kps[5][1]
    right_neck_y = template_kps[5][0]

    # print((left_neck_y,left_neck_x))
    # cv2.circle(test_img,(left_neck_x,left_neck_y),6,(243,56,32),-1)
    # cv2.circle(test_img,(right_neck_x,right_neck_y),6,(243,56,32),-1)
    # cv2.imwrite('neck.jpg', test_img)

    orig_left_neck_x,orig_left_neck_y = new_coordinates_after_resize_img((h,w), (width,height), (left_neck_x, left_neck_y)) # just modify this line
    orig_right_neck_x,orig_right_neck_y = new_coordinates_after_resize_img((h,w), (width,height), (right_neck_x, right_neck_y)) # just modify this line

    orig_w = abs(orig_left_neck_x-orig_right_neck_x)

    cv2.circle(input,(orig_left_neck_x,orig_left_neck_y),10,(243,56,32),-1)
    cv2.circle(input,(orig_right_neck_x,orig_right_neck_y),10,(243,56,32),-1)
    cv2.imwrite('orig.jpg', input)

    from mtcnn.mtcnn import MTCNN
    detector = MTCNN()
    rects = detector.detect_faces(input)
    eye = rects[0]['keypoints']['left_eye']
    dist = min(abs(orig_left_neck_y-eye[1]),orig_left_neck_y)

    y = orig_left_neck_y
    x = max(orig_left_neck_x-40,40)

    h = min(dist+int(dist//1.2),orig_left_neck_y)
    print(h,"=============")
    neck_dist = min(orig_right_neck_x+40,orig_right_neck_x+abs(orig_right_neck_x+width))
    w = abs(neck_dist-x)

    crop_img = input[y-h+20:y, x:x+w]
    cv2.imwrite("cropped.jpg", crop_img)
    return crop_img

if __name__=="__main__":
  template_path = "12.jpg"
  input = cv2.imread(template_path)
  poseNet(input)

max_issues_repo_name: aminul-palash/Background-Remove-API

max_forks_repo_name: aminul-palash/Background-Remove-API


text:
#!/usr/bin/python
# -*- coding: utf-8 -*-

"""
    Copyright 2010 Naranjo, Manuel Francisco <manuel@aircable.net>

    Licensed under the Apache License, Version 2.0 (the "License");
    you may not use this file except in compliance with the License.
    You may obtain a copy of the License at

	http://www.apache.org/licenses/LICENSE-2.0

    Unless required by applicable law or agreed to in writing, software
    distributed under the License is distributed on an "AS IS" BASIS,
    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    See the License for the specific language governing permissions and
    limitations under the License.

Based on one of the matplot examples, this will draw received data from the
serial port on real time.
"""

import os, sys, random, select, termios, fcntl, serial

from matplotlib.figure import Figure
from matplotlib.backends.backend_qt4agg import FigureCanvasQTAgg as FigureCanvas

from PyQt4 import QtCore, QtGui

import numpy as np
import time

SPAN=60 # last X seconds

def scan():
    """scan for available ports. return a list of tuples (num, name)"""
    available = []
    for i in range(256):
        try:
            s = serial.Serial(i)
            available.append((i, s.portstr))
            s.close()   # explicit close 'cause of delayed GC in java
        except serial.SerialException:
            pass
    if os.name == 'posix':
      for i in range(256):
        try:
            s = serial.Serial('/dev/ttyUSB%s' % i)
            available.append( (i, s.portstr) )
            s.close()   # explicit close 'cause of delayed GC in java
        except serial.SerialException:
            pass
    return available

class BlitQT(FigureCanvas):

    def addSerie(self, serie):
        ind = len(self.series)+1
        plot = self.figure.add_subplot( ind, 1, ind )
        plot.autoscale(enable=True, tight=True)
        plot.set_ylabel(serie)
        plot.set_xlim(left=self.tstart)
        plot.grid()

        self.series[serie] = {
    	    'plot': plot,
	    'x':[],
	    'y':[],
	    'background': self.copy_from_bbox(plot.bbox),
	    'chart': plot.plot([],[])[0]
	}

	i=1
	for s in self.series:
	    if s == serie:
		continue
	    self.figure.delaxes(self.series[s]['plot'])
	    plot = self.figure.add_subplot( ind, 1, i)
    	    plot.autoscale(enable=True, tight=True)
    	    plot.set_ylabel(s)
    	    plot.grid()
    	    self.series[s]['plot']=plot
	    self.series[s]['background'] = self.copy_from_bbox(plot.bbox)
	    self.series[s]['chart']=plot.plot([],[])[0]
	    i+=1
	self.draw()

    def addDataPoint(self, serie, value):
	if serie not in self.series:
	    # no conozco serie, agregar
	    self.addSerie(serie)

	plot = self.series[serie]
	self.lasttime = time.time()
	plot['x'].append(self.lasttime)
	plot['y'].append(value)
	plot['chart'].set_xdata(plot['x'])
	plot['chart'].set_ydata(plot['y'])
	if len(plot['x']) == 0:
	    return
	m = min(plot['y'])*0.9
	M = max(plot['y'])*1.1
	if m != M:
	    plot['plot'].set_ylim(m, M)

    def updateSpan(self):
	now = time.time()
	for s in self.series:
	    self.series[s]['plot'].set_xlim(now-SPAN, now)

    def keyPressEvent(self, event):
	key = event.key()
	if key == QtCore.Qt.Key_Escape:
	    self.close
	    return
	print key, chr(key)
	self.port.write('%s' % chr(key))

    def __init__(self, port):
        FigureCanvas.__init__(self, Figure())
        self.series = dict()
        self.draw()
        self.tstart = time.time()
        self.startTimer(500)
        self.port = serial.Serial(port, timeout=0.01)
        self.lasttime = 0
        self.cnt = 0

	# vaciar buffer
	if len(select.select([self.port,],[],[], 0)[0]) > 0:
	    self.port.flush()

    def timerEvent(self, event):
	self.cnt+=1

	# ver si hay algo en el puerto
	if len(select.select([self.port,],[],[], 0)[0]) == 0:
	    self.doRedraw()
	    return

	for line in self.port.readlines():
	    line=line.strip()
	    if len(line) == 0:
		continue
	    try:
		serie, dato = line.split(';')[:2]
		dato = float(dato)
		self.addDataPoint(serie, dato)
	    except Exception, err:
		pass
	self.doRedraw()

    def resizeEvent(self, event):
	for i in self.series:
	    self.series[i]['background']=self.copy_from_bbox(self.series[i]['plot'].bbox)
	super(BlitQT, self).resizeEvent(event)

    def doRedraw(self):
	self.updateSpan()
	self.draw()
	for i in self.series:
	    pl = self.series[i]
	    # clear!
	    self.restore_region(['background'])
	    # update chart
	    pl['plot'].draw_artist(pl['chart'])
	    # plot axis
	    self.blit(pl['plot'].bbox)
	self.draw()

if __name__=='__main__':
    i = 0
    ports = scan()
    if len(ports) == 0 and len(sys.argv) == 1:
	print "Couldn't found port, please give me a port name as argument"
	sys.exit(0)

    if len(sys.argv) == 1:
	print "Choose port:"
	for num, port in ports:
	    print "\t %i: %s" % (i, port)
	    i+=1
	port = int(raw_input("enter number: "))
	port = ports[port][1]
    else:
	port = sys.argv[1]

    app = QtGui.QApplication(sys.argv)
    widget = BlitQT(port)
    widget.show()

    sys.exit(app.exec_())

max_stars_repo_name: manuelnaranjo/pyrealtimecharter

max_issues_repo_name: manuelnaranjo/pyrealtimecharter

max_forks_repo_name: manuelnaranjo/pyrealtimecharter


text:
import sys
sys.path.append('../')

from pathlib import Path
import pickle
import matplotlib.pyplot as plt
import numpy as np

from py_diff_pd.common.common import print_info, print_error

from py_diff_pd.common.common import ndarray

if __name__ == '__main__':
    folder = Path('routing_tendon_3d')
    import multiprocessing
    cpu_cnt = multiprocessing.cpu_count()
    thread_ct = cpu_cnt - 1
    try:
        data = pickle.load(open(folder / 'data_{:04d}_threads.bin'.format(thread_ct), 'rb'))
    except:
        print_error('Log file not found.')
    loss_l, loss_h = data['loss_range']
    print_info('Loss range: {:3f}, {:3f}'.format(loss_l, loss_h))
    def normalize_loss(unnormalized_loss):
        return (unnormalized_loss - loss_l) / (loss_h - loss_l)

    data_file = Path('routing_tendon_3d') / 'data_{:04d}_threads.bin'.format(thread_ct)
    if data_file.exists():
        print_info('Loading {}'.format(data_file))
        data = pickle.load(open(data_file, 'rb'))
        for method in ['newton_pcg', 'newton_cholesky', 'pd_eigen']:
            total_time = 0
            avg_forward = 0
            average_backward = 0
            for d in data[method]:
                d['loss'] = normalize_loss(d['loss'])
                print('loss: {:8.3f}, |grad|: {:8.3f}, forward time: {:6.3f}s, backward time: {:6.3f}s'.format(
                    d['loss'], np.linalg.norm(d['grad']), d['forward_time'], d['backward_time']))
                total_time += d['forward_time'] + d['backward_time']
                average_backward += d['backward_time']
                avg_forward += d['forward_time']
            avg_forward /= len(data[method])
            average_backward /= len(data[method])
            print_info('Optimizing with {} finished in {:6.3f}s with {:d} iterations. Average Backward time: {:6.3f}s, Average Forward Time = {:6.3f}s'.format(
                method, total_time,  len(data[method]), average_backward, avg_forward))

    plt.rc('pdf', fonttype=42)
    plt.rc('font', size=30)             # Controls default text sizes.
    plt.rc('axes', titlesize=36)        # Fontsize of the axes title.
    plt.rc('axes', labelsize=36)        # Fontsize of the x and y labels.
    plt.rc('xtick', labelsize=36)       # Fontsize of the tick labels.
    plt.rc('ytick', labelsize=36)       # Fontsize of the tick labels.
    plt.rc('legend', fontsize=36)       # Legend fontsize.
    plt.rc('figure', titlesize=36)      # Fontsize of the figure title.

    acts = {}
    losses = {}
    for method in ['newton_pcg', 'newton_cholesky', 'pd_eigen']:
        acts[method] = [np.linalg.norm(d['x']) for d in data[method]]
        losses[method] = [d['loss'] for d in data[method]]

    fig = plt.figure(figsize=(20, 10))

    ax_act = fig.add_subplot(121)

    ax_loss= fig.add_subplot(122)

    titles = ['muscle actuation', 'loss']
    for title, ax, y in zip(titles, (ax_act, ax_loss), (acts, losses)):

        if 'muscle' in title:
            ax.set_ylabel("|actuation|")
            ax.grid(True, which='both')
        else:
            ax.set_ylabel("loss")
            ax.set_yscale('log')
            ax.grid(True)
        ax.set_xlabel('function evaluations')
        for method, method_ref_name, color in zip(['newton_pcg', 'newton_cholesky', 'pd_eigen'],
            ['PCG', 'Cholesky', 'Ours'], ['tab:blue', 'tab:red', 'tab:green']):
            ax.plot(y[method], color=color, label=method_ref_name, linewidth=4)
        ax.set_title(title, pad=25)
        ax.set_xticks([0, 10, 20, 30, 40])
        handles, labels = ax.get_legend_handles_labels()

    plt.subplots_adjust(bottom = 0.25, wspace=0.3)
    # Share legends.
    fig.legend(handles, labels, loc='lower center', ncol=3)#, bbox_to_anchor=(0.5, 0.17))

    fig.savefig(folder / 'tendon.pdf')
    plt.show()


text:
"""
Copyright 2020 The OneFlow Authors. All rights reserved.

Licensed under the Apache License, Version 2.0 (the "License");
you may not use this file except in compliance with the License.
You may obtain a copy of the License at

    http://www.apache.org/licenses/LICENSE-2.0

Unless required by applicable law or agreed to in writing, software
distributed under the License is distributed on an "AS IS" BASIS,
WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
See the License for the specific language governing permissions and
limitations under the License.
"""
import unittest
import os
import numpy as np
import tempfile

import oneflow as flow
import oneflow.unittest

def _test_linear_graph_save_load(test_case, device):
    def train_with_graph(call_cnt=0, state_dict_dir=None, last_state_dict=None):
        linear = flow.nn.Linear(3, 8)
        linear = linear.to(device)
        flow.nn.init.constant_(linear.weight, 2.068758)
        flow.nn.init.constant_(linear.bias, 0.23)
        of_sgd = flow.optim.SGD(linear.parameters(), lr=0.001, momentum=0.9)

        x = flow.tensor(
            [
                [-0.94630778, -0.83378579, -0.87060891],
                [2.0289922, -0.28708987, -2.18369248],
                [0.35217619, -0.67095644, -1.58943879],
                [0.08086036, -1.81075924, 1.20752494],
                [0.8901075, -0.49976737, -1.07153746],
                [-0.44872912, -1.07275683, 0.06256855],
                [-0.22556897, 0.74798368, 0.90416439],
                [0.48339456, -2.32742195, -0.59321527],
            ],
            dtype=flow.float32,
            device=device,
            requires_grad=False,
        )

        class LinearTrainGraph(flow.nn.Graph):
            def __init__(self):
                super().__init__()
                self.linear = linear
                self.add_optimizer(of_sgd)

            def build(self, x):
                out = self.linear(x)
                out = out.sum()
                out.backward()
                return out

        linear_t_g = LinearTrainGraph()
        if call_cnt == 1:
            state_dict = flow.load(state_dict_dir)
            linear_t_g.load_state_dict(state_dict)
            # Check state in module has been loaded.
            test_case.assertTrue(
                np.array_equal(state_dict["linear"]["weight"].numpy(), linear.weight)
            )
            test_case.assertTrue(
                np.array_equal(state_dict["linear"]["bias"].numpy(), linear.bias)
            )
        # Get state dict before compile is allowed.
        init_state_dict = linear_t_g.state_dict()

        of_graph_out = linear_t_g(x)
        iter0_state_dict = linear_t_g.state_dict()
        if call_cnt == 1:
            # Check additional variable state initialized in job has been loaded.
            cur_train_step = iter0_state_dict["System-Train-TrainStep"].numpy()[0]
            test_case.assertEqual(3, cur_train_step)
            test_case.assertTrue(
                cur_train_step == last_state_dict["System-Train-TrainStep"].numpy()[0]
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear"]["weight"].numpy(),
                    last_state_dict["linear"]["weight"].numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear"]["bias"].numpy(),
                    last_state_dict["linear"]["bias"].numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear.weight-momentum"].numpy(),
                    last_state_dict["linear.weight-momentum"].numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear.bias-momentum"].numpy(),
                    last_state_dict["linear.bias-momentum"].numpy(),
                )
            )

        of_graph_out = linear_t_g(x)
        of_graph_out.numpy()
        iter1_state_dict = linear_t_g.state_dict()
        if call_cnt == 0:
            flow.save(iter1_state_dict, state_dict_dir)

        if call_cnt == 0:
            of_graph_out = linear_t_g(x)
            iter2_state_dict = linear_t_g.state_dict()
            of_graph_out.numpy()
            return iter2_state_dict

    with tempfile.TemporaryDirectory(prefix="graph_save_load_local") as state_dict_dir:
        iter2_state_dict = train_with_graph(0, state_dict_dir)
        train_with_graph(1, state_dict_dir, iter2_state_dict)

@unittest.skipIf(os.getenv("ONEFLOW_TEST_CPU_ONLY"), "only test cpu cases")
@flow.unittest.skip_unless_1n1d()
class TestLinearGraphSaveLoad(oneflow.unittest.TestCase):
    def test_linear_graph_save_load_gpu(test_case):
        _test_linear_graph_save_load(test_case, flow.device("cuda"))

    def _test_linear_graph_save_load_cpu(test_case):
        _test_linear_graph_save_load(test_case, flow.device("cpu"))

def _test_linear_graph_save_load_global(test_case, device):
    P = flow.placement("cuda", ranks=[0, 1])
    B = flow.sbp.broadcast
    S = flow.sbp.split(0)

    def train_with_graph(call_cnt=0, state_dict_dir=None, last_state_dict=None):
        linear = flow.nn.Linear(3, 8)
        linear = linear.to(device)
        flow.nn.init.constant_(linear.weight, 2.068758)
        flow.nn.init.constant_(linear.bias, 0.23)
        linear.to_global(placement=P, sbp=B)
        of_sgd = flow.optim.SGD(linear.parameters(), lr=0.001, momentum=0.9)

        x = flow.tensor(
            [
                [-0.94630778, -0.83378579, -0.87060891],
                [2.0289922, -0.28708987, -2.18369248],
                [0.35217619, -0.67095644, -1.58943879],
                [0.08086036, -1.81075924, 1.20752494],
                [0.8901075, -0.49976737, -1.07153746],
                [-0.44872912, -1.07275683, 0.06256855],
                [-0.22556897, 0.74798368, 0.90416439],
                [0.48339456, -2.32742195, -0.59321527],
            ],
            dtype=flow.float32,
            device=device,
            requires_grad=False,
        )
        x = x.to_global(placement=P, sbp=S)

        class LinearTrainGraphGlobal(flow.nn.Graph):
            def __init__(self):
                super().__init__()
                self.linear = linear
                self.add_optimizer(of_sgd)

            def build(self, x):
                out = self.linear(x)
                out = out.sum()
                out.backward()
                return out

        linear_t_g = LinearTrainGraphGlobal()
        if call_cnt == 1:
            state_dict = flow.load(state_dict_dir, global_src_rank=0)
            linear_t_g.load_state_dict(state_dict)
            # Check state in module has been loaded.
            # Tensors in state dict are save to rank 0, so they need to be broadcast to rank 0 and 1 before check.
            test_case.assertTrue(
                np.array_equal(
                    state_dict["linear"]["weight"]
                    .to_global(placement=P, sbp=B)
                    .to_local()
                    .numpy(),
                    linear.weight.to_local().numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    state_dict["linear"]["bias"]
                    .to_global(placement=P, sbp=B)
                    .to_local()
                    .numpy(),
                    linear.bias.to_local().numpy(),
                )
            )
        # Get state dict before compile is allowed.
        init_state_dict = linear_t_g.state_dict()

        of_graph_out = linear_t_g(x)
        iter0_state_dict = linear_t_g.state_dict()
        if call_cnt == 1:
            # Check additional variable state initialized in job has been loaded.
            # TrainStep's placement is only on rank 0, so it needs to be broadcast to rank 0 and 1 before check.
            cur_train_step = (
                iter0_state_dict["System-Train-TrainStep"]
                .to_global(placement=P, sbp=B)
                .to_local()
                .numpy()[0]
            )
            test_case.assertEqual(3, cur_train_step)
            test_case.assertTrue(
                cur_train_step
                == last_state_dict["System-Train-TrainStep"]
                .to_global(placement=P, sbp=B)
                .to_local()
                .numpy()[0]
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear"]["weight"].to_local().numpy(),
                    last_state_dict["linear"]["weight"].to_local().numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear"]["bias"].to_local().numpy(),
                    last_state_dict["linear"]["bias"].to_local().numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear.weight-momentum"].to_local().numpy(),
                    last_state_dict["linear.weight-momentum"].to_local().numpy(),
                )
            )
            test_case.assertTrue(
                np.array_equal(
                    iter0_state_dict["linear.bias-momentum"].to_local().numpy(),
                    last_state_dict["linear.bias-momentum"].to_local().numpy(),
                )
            )

        of_graph_out = linear_t_g(x)
        of_graph_out.numpy()
        iter1_state_dict = linear_t_g.state_dict()
        if call_cnt == 0:
            flow.save(iter1_state_dict, state_dict_dir, global_dst_rank=0)

        if call_cnt == 0:
            of_graph_out = linear_t_g(x)
            of_graph_out.numpy()
            iter2_state_dict = linear_t_g.state_dict()
            return iter2_state_dict

    with tempfile.TemporaryDirectory(prefix="graph_save_load_global") as state_dict_dir:
        iter2_state_dict = train_with_graph(0, state_dict_dir)
        train_with_graph(1, state_dict_dir, iter2_state_dict)

@unittest.skipIf(os.getenv("ONEFLOW_TEST_CPU_ONLY"), "only test cpu cases")
@flow.unittest.skip_unless_1n2d()
class TestLinearGraphSaveLoadGlobal(oneflow.unittest.TestCase):
    def test_linear_graph_save_load_gpu(test_case):
        _test_linear_graph_save_load_global(test_case, flow.device("cuda"))

    def _test_linear_graph_save_load_cpu(test_case):
        _test_linear_graph_save_load_global(test_case, flow.device("cpu"))

if __name__ == "__main__":
    unittest.main()


text:
# rFirm
# See full license in LICENSE.txt.

import os
import logging

import numpy as np
import pandas as pd

from activitysim.core import inject
from activitysim.core.config import setting

logger = logging.getLogger(__name__)

def scenario_dir():

    scenarios_dir = setting('scenarios_dir')
    assert scenarios_dir is not None, "scenarios_dir not defined in settings file"

    if not os.path.exists(scenarios_dir):
        raise RuntimeError("scenarios_dir not found: %s" % scenarios_dir)

    scenario_name = setting('scenario_name')
    assert scenario_name is not None, "scenario_name not defined in settings file"

    scenario_dir_path = os.path.join(scenarios_dir, scenario_name)
    assert os.path.exists(scenario_dir_path), "scenario_dir not found: %s" % scenario_dir_path

    return scenario_dir_path

def read_table(table_name, table_info, data_dir):

    logger.info("read_table read %s" % table_name)

    # read the csv file
    data_filename = table_info.get('filename', '%s.csv' % table_name)
    data_file_path = os.path.join(data_dir, data_filename)
    if not os.path.exists(data_file_path):
        raise RuntimeError("read_table %s - input file not found: %s"
                           % (table_name, data_file_path,))

    logger.info("Reading csv file %s" % data_file_path)
    df = pd.read_csv(data_file_path, comment='#')

    # rename columns
    column_map = table_info.get('column_map', None)
    if column_map:
        df.rename(columns=column_map, inplace=True)

    # set index
    index_col = table_info.get('index_col', None)
    if index_col is not None:
        if index_col in df.columns:
            df.set_index(index_col, inplace=True)
        else:
            df.index.names = [index_col]

    return df

def load_tables(table_list_name, data_dir=None):

    table_list = setting(table_list_name)
    if table_list is None:
        raise "I expected to find table list '%s' with table_info in settings." % table_list_name

    if data_dir is None:
        data_dir = setting('data_dir', inject.get_injectable('data_dir'))

    for table_name, table_info in table_list.iteritems():

        df = read_table(table_name, table_info, data_dir)
        inject.add_table(table_name, df)

def bucket_round(x, threshold=0.5):

    n = len(x)
    vecf = x.astype(int)
    vecd = x % 1.0
    veca = [0] * n
    adj = 0
    for i in range(len(x)):
        adj += vecd[i]
        if adj >= threshold:
            veca[i] = 1
            adj -= 1

    return (vecf + veca)

def target_round(weights, target_sum=None):
    """
    Round weights while ensuring (as far as possible that result sums to target_sum)

    Parameters
    ----------
    weights : numpy.ndarray(float)
    target_sum : int

    Returns
    -------
    rounded_weights : numpy.ndarray array of ints
    """

    if target_sum is None:
        target_sum = round(weights.sum())

    assert target_sum == int(target_sum)

    target_sum = int(target_sum)

    # integer part of numbers to round (astype both copies and coerces)
    rounded_weights = weights.astype(int)
    resid_weights = weights % 1.0

    # find number of residuals that we need to round up
    int_shortfall = target_sum - rounded_weights.sum()

    # clip to feasible, in case target was not achievable by rounding
    int_shortfall = np.clip(int_shortfall, 0, len(resid_weights))

    # Order the residual weights and round at the tipping point where target_sum is achieved
    if int_shortfall > 0:
        # indices of the int_shortfall highest resid_weights
        i = np.argsort(resid_weights)[-int_shortfall:]

        # add 1 to the integer weights that we want to round upwards
        rounded_weights[i] += 1

    assert rounded_weights.sum() == target_sum

    return rounded_weights

max_stars_repo_name: RSGInc/rFirm

max_issues_repo_name: RSGInc/rFirm

max_forks_repo_name: RSGInc/rFirm


text:
import numpy as np

#pt except for reshaping, then discriminant
btag_feval_dims = {0:[1],1:[1],2:[1],3:[2]}

def convert_btag_csv_file(csvFilePath):
    btag_f = open(csvFilePath)
    nameandcols = btag_f.readline().split(';')
    btag_f.close()
    name = nameandcols[0].strip()
    columns = nameandcols[1].strip()
    columns = [column.strip() for column in columns.split(',')]

    corrections = np.genfromtxt(csvFilePath,
                                dtype=None,
                                names=tuple(columns),
                                converters={1:lambda s: s.strip(),
                                            2:lambda s: s.strip(),
                                           10:lambda s: s.strip(' "')},
                                delimiter = ',',
                                skip_header=1,
                                unpack=True,
                                encoding='ascii'
                                )

    all_names = corrections[[columns[i] for i in range(4)]]
    labels = np.unique(corrections[[columns[i] for i in range(4)]])
    names_and_bins = np.unique(corrections[[columns[i] for i in [0,1,2,3,4,6,8]]])
    wrapped_up = {}
    for label in labels:
        etaMins = np.unique(corrections[np.where(all_names == label)][columns[4]])
        etaMaxs = np.unique(corrections[np.where(all_names == label)][columns[5]])
        etaBins = np.union1d(etaMins,etaMaxs)
        ptMins = np.unique(corrections[np.where(all_names == label)][columns[6]])
        ptMaxs = np.unique(corrections[np.where(all_names == label)][columns[7]])
        ptBins = np.union1d(ptMins,ptMaxs)
        discrMins = np.unique(corrections[np.where(all_names == label)][columns[8]])
        discrMaxs = np.unique(corrections[np.where(all_names == label)][columns[9]])
        discrBins = np.union1d(discrMins,discrMaxs)
        vals = np.zeros(shape=(len(discrBins)-1,len(ptBins)-1,len(etaBins)-1),
                        dtype=corrections.dtype[10])
        for i,eta_bin in enumerate(etaBins[:-1]):
            for j,pt_bin in enumerate(ptBins[:-1]):
                for k,discr_bin in enumerate(discrBins[:-1]):
                    this_bin = np.where((all_names == label) &
                                        (corrections[columns[4]] == eta_bin) &
                                        (corrections[columns[6]] == pt_bin)  &
                                        (corrections[columns[8]] == discr_bin))
                    vals[k,j,i] = corrections[this_bin][columns[10]][0]
        label_decode = []
        for i in range(len(label)):
            label_decode.append(label[i])
            if isinstance(label_decode[i],bytes):
                label_decode[i] = label_decode[i].decode()
            else:
                label_decode[i] = str(label_decode[i])
        str_label = '_'.join([name]+label_decode)
        feval_dim = btag_feval_dims[label[0]]
        wrapped_up[(str_label,'dense_evaluated_lookup')] = (vals,(etaBins,ptBins,discrBins),tuple(feval_dim))
    return wrapped_up


text:
##
# \file reconstruct_volume_from_slices.py
# \brief      Script to reconstruct an isotropic, high-resolution volume from
#             multiple motion-corrected (or static) stacks of low-resolution 2D
#             slices.
#
# \author     Michael Ebner (michael.ebner.14@ucl.ac.uk)
# \date       March 2017
#

import os
import numpy as np
import SimpleITK as sitk

import pysitk.python_helper as ph
import pysitk.simple_itk_helper as sitkh

import niftymic.base.stack as st
import niftymic.base.data_reader as dr
import niftymic.base.data_writer as dw
import niftymic.reconstruction.admm_solver as admm
import niftymic.utilities.intensity_correction as ic
import niftymic.reconstruction.primal_dual_solver as pd
import niftymic.reconstruction.tikhonov_solver as tk
import niftymic.reconstruction.scattered_data_approximation as sda
import niftymic.utilities.binary_mask_from_mask_srr_estimator as bm
from niftymic.utilities.input_arparser import InputArgparser

from niftymic.definitions import ALLOWED_EXTENSIONS

def main():

    time_start = ph.start_timing()

    # Set print options for numpy
    np.set_printoptions(precision=3)

    # Read input
    input_parser = InputArgparser(
        description="Volumetric MRI reconstruction framework to reconstruct "
        "an isotropic, high-resolution 3D volume from multiple "
        "motion-corrected (or static) stacks of low-resolution slices.",
    )
    input_parser.add_filenames(required=True)
    input_parser.add_filenames_masks()
    input_parser.add_dir_input_mc()
    input_parser.add_output(required=True)
    input_parser.add_suffix_mask(default="_mask")
    input_parser.add_target_stack(default=None)
    input_parser.add_extra_frame_target(default=10)
    input_parser.add_isotropic_resolution(default=None)
    input_parser.add_intensity_correction(default=1)
    input_parser.add_reconstruction_space(default=None)
    input_parser.add_minimizer(default="lsmr")
    input_parser.add_iter_max(default=10)
    input_parser.add_reconstruction_type(default="TK1L2")
    input_parser.add_data_loss(default="linear")
    input_parser.add_data_loss_scale(default=1)
    input_parser.add_alpha(
        default=0.01  # TK1L2 @ isotropic_resolution = 0.8
        # default=0.006  #TVL2, HuberL2 @ isotropic_resolution = 0.8
    )
    input_parser.add_rho(default=0.1)
    input_parser.add_tv_solver(default="PD")
    input_parser.add_pd_alg_type(default="ALG2")
    input_parser.add_iterations(default=15)
    input_parser.add_log_config(default=1)
    input_parser.add_use_masks_srr(default=0)
    input_parser.add_slice_thicknesses(default=None)
    input_parser.add_verbose(default=0)
    input_parser.add_viewer(default="itksnap")
    input_parser.add_argument(
        "--mask", "-mask",
        action='store_true',
        help="If given, input images are interpreted as image masks. "
        "Obtained volumetric reconstruction will be exported in uint8 format."
    )
    input_parser.add_argument(
        "--sda", "-sda",
        action='store_true',
        help="If given, the volume is reconstructed using "
        "Scattered Data Approximation (Vercauteren et al., 2006). "
        "--alpha is considered the value for the standard deviation then. "
        "Recommended value is, e.g., --alpha 0.8"
    )

    args = input_parser.parse_args()
    input_parser.print_arguments(args)

    if args.reconstruction_type not in ["TK1L2", "TVL2", "HuberL2"]:
        raise IOError("Reconstruction type unknown")

    if np.alltrue([not args.output.endswith(t) for t in ALLOWED_EXTENSIONS]):
        raise ValueError(
            "output filename '%s' invalid; "
            "allowed image extensions are: %s" % (
                args.output, ", ".join(ALLOWED_EXTENSIONS)))

    dir_output = os.path.dirname(args.output)
    ph.create_directory(dir_output)

    debug = 0

    if args.log_config:
        input_parser.log_config(os.path.abspath(__file__))

    if args.verbose:
        show_niftis = []
        # show_niftis = [f for f in args.filenames]

    # --------------------------------Read Data--------------------------------
    ph.print_title("Read Data")

    if args.mask:
        filenames_masks = args.filenames
    else:
        filenames_masks = args.filenames_masks

    data_reader = dr.MultipleImagesReader(
        file_paths=args.filenames,
        file_paths_masks=filenames_masks,
        suffix_mask=args.suffix_mask,
        dir_motion_correction=args.dir_input_mc,
        stacks_slice_thicknesses=args.slice_thicknesses,
    )
    data_reader.read_data()
    stacks = data_reader.get_data()

    ph.print_info("%d input stacks read for further processing" % len(stacks))

    # Specify target stack for intensity correction and reconstruction space
    if args.target_stack is None:
        target_stack_index = 0
    else:
        # TODO: deal with case when target stack got rejected in previous step
        filenames = ["%s.nii.gz" % s.get_filename() for s in stacks]
        filename_target_stack = os.path.basename(args.target_stack)
        try:
            target_stack_index = filenames.index(filename_target_stack)
        except ValueError as e:
            raise ValueError(
                "--target-stack must correspond to an image as provided by "
                "--filenames")

    # ---------------------------Intensity Correction--------------------------
    if args.intensity_correction and not args.mask:
        ph.print_title("Intensity Correction")
        intensity_corrector = ic.IntensityCorrection()
        intensity_corrector.use_individual_slice_correction(False)
        intensity_corrector.use_stack_mask(True)
        intensity_corrector.use_reference_mask(True)
        intensity_corrector.use_verbose(False)

        for i, stack in enumerate(stacks):
            if i == target_stack_index:
                ph.print_info("Stack %d (%s): Reference image. Skipped." % (
                    i + 1, stack.get_filename()))
                continue
            else:
                ph.print_info("Stack %d (%s): Intensity Correction ... " % (
                    i + 1, stack.get_filename()), newline=False)
            intensity_corrector.set_stack(stack)
            intensity_corrector.set_reference(
                stacks[target_stack_index].get_resampled_stack(
                    resampling_grid=stack.sitk,
                    interpolator="NearestNeighbor",
                ))
            intensity_corrector.run_linear_intensity_correction()
            stacks[i] = intensity_corrector.get_intensity_corrected_stack()
            print("done (c1 = %g) " %
                  intensity_corrector.get_intensity_correction_coefficients())

    # -------------------------Volumetric Reconstruction-----------------------
    ph.print_title("Volumetric Reconstruction")

    # Reconstruction space defined by isotropically resampled,
    # bounding box-cropped target stack
    if args.reconstruction_space is None:
        recon0 = stacks[target_stack_index].get_isotropically_resampled_stack(
            resolution=args.isotropic_resolution,
            extra_frame=args.extra_frame_target,
        )
        recon0 = recon0.get_cropped_stack_based_on_mask(
            boundary_i=args.extra_frame_target,
            boundary_j=args.extra_frame_target,
            boundary_k=args.extra_frame_target,
            unit="mm",
        )

    # Reconstruction space was provided by user
    else:
        recon0 = st.Stack.from_filename(args.reconstruction_space,
                                        extract_slices=False)

        # Change resolution for isotropic resolution if provided by user
        if args.isotropic_resolution is not None:
            recon0 = recon0.get_isotropically_resampled_stack(
                args.isotropic_resolution)

        # Use image information of selected target stack as recon0 serves
        # as initial value for reconstruction
        recon0 = stacks[target_stack_index].get_resampled_stack(recon0.sitk)
        recon0 = recon0.get_stack_multiplied_with_mask()

    ph.print_info(
        "Reconstruction space defined with %s mm3 resolution" %
        " x ".join(["%.2f" % s for s in recon0.sitk.GetSpacing()])
    )

    if debug:
        # visualize (intensity corrected) data alongside recon0 init
        show = [st.Stack.from_stack(s) for s in stacks]
        show.insert(0, recon0)
        sitkh.show_stacks(show)

    if args.sda:
        ph.print_title("Compute SDA reconstruction")
        SDA = sda.ScatteredDataApproximation(
            stacks, recon0, sigma=args.alpha, sda_mask=args.mask)
        SDA.run()
        recon = SDA.get_reconstruction()
        filename = SDA.get_setting_specific_filename()
        if args.mask:
            dw.DataWriter.write_mask(
                recon.sitk_mask, args.output, description=filename)
        else:
            dw.DataWriter.write_image(
                recon.sitk, args.output, description=filename)

        if args.verbose:
            show_niftis.insert(0, args.output)

    else:
        if args.reconstruction_type in ["TVL2", "HuberL2"]:
            ph.print_title(
                "Compute Initial value for %s" % args.reconstruction_type)
            SRR0 = sda.ScatteredDataApproximation(stacks, recon0, sigma=0.8)
        else:
            ph.print_title(
                "Compute %s reconstruction" % args.reconstruction_type)
            SRR0 = tk.TikhonovSolver(
                stacks=stacks,
                reconstruction=recon0,
                alpha=args.alpha,
                iter_max=args.iter_max,
                reg_type="TK1",
                minimizer=args.minimizer,
                data_loss=args.data_loss,
                data_loss_scale=args.data_loss_scale,
                use_masks=args.use_masks_srr,
                # verbose=args.verbose,
            )
        SRR0.run()

        recon = SRR0.get_reconstruction()
        filename = SRR0.get_setting_specific_filename()

        if args.verbose and args.reconstruction_type in ["TVL2", "HuberL2"]:
            output = ph.append_to_filename(args.output, "_init")

            if args.mask:
                mask_estimator = bm.BinaryMaskFromMaskSRREstimator(recon.sitk)
                mask_estimator.run()
                mask_sitk = mask_estimator.get_mask_sitk()
                dw.DataWriter.write_mask(
                    mask_sitk, output, description=filename)
            else:
                dw.DataWriter.write_image(
                    recon.sitk, output, description=filename)

            show_niftis.insert(0, output)

        if args.reconstruction_type in ["TVL2", "HuberL2"]:
            ph.print_title("Compute %s reconstruction" %
                           args.reconstruction_type)
            if args.tv_solver == "ADMM":
                SRR = admm.ADMMSolver(
                    stacks=stacks,
                    reconstruction=st.Stack.from_stack(
                        SRR0.get_reconstruction()),
                    minimizer=args.minimizer,
                    alpha=args.alpha,
                    iter_max=args.iter_max,
                    rho=args.rho,
                    data_loss=args.data_loss,
                    iterations=args.iterations,
                    use_masks=args.use_masks_srr,
                    verbose=args.verbose,
                )

            else:
                SRR = pd.PrimalDualSolver(
                    stacks=stacks,
                    reconstruction=st.Stack.from_stack(
                        SRR0.get_reconstruction()),
                    minimizer=args.minimizer,
                    alpha=args.alpha,
                    iter_max=args.iter_max,
                    iterations=args.iterations,
                    alg_type=args.pd_alg_type,
                    reg_type="TV" if args.reconstruction_type == "TVL2" else "huber",
                    data_loss=args.data_loss,
                    use_masks=args.use_masks_srr,
                    verbose=args.verbose,
                )
            SRR.run()
            recon = SRR.get_reconstruction()
            filename = SRR.get_setting_specific_filename()

        if args.mask:
            mask_estimator = bm.BinaryMaskFromMaskSRREstimator(recon.sitk)
            mask_estimator.run()
            mask_sitk = mask_estimator.get_mask_sitk()
            dw.DataWriter.write_mask(
                mask_sitk, args.output, description=filename)

        else:
            dw.DataWriter.write_image(
                recon.sitk, args.output, description=filename)

        if args.verbose:
            show_niftis.insert(0, args.output)

    if args.verbose:
        ph.show_niftis(show_niftis, viewer=args.viewer)

    ph.print_line_separator()

    elapsed_time = ph.stop_timing(time_start)
    ph.print_title("Summary")
    exe_file_info = os.path.basename(os.path.abspath(__file__)).split(".")[0]
    print("%s | Computational Time for Volumetric Reconstruction: %s" % (
        exe_file_info, elapsed_time))

    return 0

if __name__ == '__main__':
    main()

max_stars_repo_name: martaranzini/NiftyMIC

max_issues_repo_name: martaranzini/NiftyMIC

max_forks_repo_name: martaranzini/NiftyMIC


text:
import matplotlib.pyplot as plt
# import pymc as pm
import numpy as np
from pertussis import a_l, a_u, J

plt.style.use('ggplot')

def draw_no_split(x, y):
    # l = y.shape[1]
    l = len(y)
    fig, axs = plt.subplots(l, figsize=(12, 8))

    return fig, axs

def draw_split(x, y, which):
    l = len(y)
    fig = plt.figure(figsize=(14, 7))
    gs = plt.GridSpec(4, 2)
    axs = []
    axs.append(plt.subplot(gs[2, 0]))
    axs.append(plt.subplot(gs[2, 1]))
    axs.append(plt.subplot(gs[0:2, 0:2]))  # Number 2 - Main Graph
    axs.append(plt.subplot(gs[3, 0]))
    axs.append(plt.subplot(gs[3, 1]))
    if which != 2:
        axs[2], axs[which] = axs[which], axs[2]
    return fig, axs

def draw_model(x, y, labels=None, split=None, collapse=False, ages = np.arange(0,J,1)):
    l = len(y)
    if l > 5:
        split = False
    if type(split) is int:
        assert l <= 5, "Must have at most 5 variables"
        fig, axs = draw_split(x, y, which=split)
    else:
        fig, axs = draw_no_split(x, y)
    for i in range(l):
        try:
            ax = axs[i]
        except:
            ax = axs
        z = y[i][ages, :]
        if collapse:
            z = y[i][:, :]
            z = z.sum(axis=0)
        lines = ax.plot(x, z.T, lw=1.3)
        try:
            ax.set_title(labels[i])
        except:
            ax.set_title("")
        ub = max(0.00001, z.max() * 1.05)
        # ub = 1
        ax.set_ylim([0, ub])
        # ax.set_xlim(left=1948)#, right = 2010)
        ax.vlines([1957, 2002, 1948], 0, ub, linestyle='--')
        if l <= 3:
            ax.legend(lines, ["{:.2f}-{:.2f}".format(a_l[j], a_u[j]) for j in ages], loc='lower left')
    if split:
        axs[split].legend(lines, [j + 1 for j in range(len(lines))], loc='lower left')

    return fig, axs

# MCMC
def plot_stoch_vars(mcmc, which=None, exclude=None):
    # Plot the variables
    # Initialize data frame for later analysis
    assert isinstance(mcmc, pm.MCMC), "mcmc should be MCMC object"
    if which:
        stoch = sorted(which)
    else:
        stoch = sorted([str(v) for v in mcmc.stochastics])
    if exclude:
        stoch.remove(exclude)
    from scipy.stats.mstats import mquantiles

    tr_len = mcmc.trace(stoch[0])[:, None].shape[0]
    # Height is determined by number of variables
    h = 3 * len(stoch)
    fig, axs = plt.subplots(len(stoch), 2, figsize=(16, h))

    for i, tr in enumerate(stoch):
        try:
            # print str(tr)
            # Prepare Values
            tr_val = mcmc.trace(tr)[:, None]
            quants = mquantiles(tr_val, prob=[0.025, 0.25, 0.5, 0.75, 0.975])

            # Plot
            ### Left: Histogram
            h = axs[i, 0].hist(tr_val, histtype='stepfilled',
                               bins=50, label=str(tr), alpha=0.6)
            m = max(h[0])
            axs[i, 0].fill_betweenx([0, m + 5], quants[0], quants[-1], color='lightgreen')
            axs[i, 0].fill_betweenx([0, m + 5], quants[1], quants[-2], color='darkgreen')
            axs[i, 0].set_title("{var} = {value:.4f} [{lci:.4f}, {hci:.4f}] ({tmp:.4f})".format(var=str(tr),
                                                                                tmp=tr_val.mean(),
                                                                                value=quants[2],
                                                                                lci=quants[0],
                                                                                hci=quants[-1]))
            axs[i, 0].set_xticks(axs[i, 0].get_xticks()[::2])  # X 0
            axs[i, 0].set_yticks(axs[i, 0].get_yticks()[::2])  # Y 0
            axs[i, 0].set_ylim([0, m + 5])
            # print axs[i,1].get_yticks()[::2]
            ### Right: Trace
            l = range(len(tr_val))
            axs[i, 1].plot(l, tr_val)
            axs[i, 1].set_xticks(axs[i, 1].get_xticks()[::2])  # X 1
            axs[i, 1].set_yticks(axs[i, 1].get_yticks()[::2])  # Y 1
            axs[i, 1].set_xlim([0, l[-1]])

        except:
            print(str(tr), " excluded")
    # fig.set_tight_layout
    fig.suptitle("Posterior Distribution <|> Convergence")
    return fig, axs

def mu_chart(mu, data):
    from pertussis import sc_ages
    from scipy.stats.mstats import mquantiles
    n_iterations, n_months, n_groups = mu.shape # Assumes (iteration, month, group)
    x = 1998 + np.arange(0, n_months, 1) / 12
    fig, axs = plt.subplots(3, int(np.ceil(n_groups/3)), figsize=(16, 16)) # Create J axes with 4 in a row
    axs = np.hstack(axs)
    for i in range(n_groups):
        ax = axs[i]
        q_mu = mu[:, :, i] # Takes group and all its lines
        quants = mquantiles(q_mu, prob=[0.05, 0.95], axis=0) #5% 95% quantiles
        ax.fill_between(x, quants[0], quants[1], color='red', alpha=0.3, label="95% CI")
        ax.plot(x, q_mu.mean(axis=0), label='Model')
        # ax.scatter(x, data[:, i], label="Data")
        ax.plot(x, data[:, i], '--', label="Data")
        ax.set_title('Age: {:.2f} - {:.2f}'.format(sc_ages[i], sc_ages[i+1]))
    axs[0].legend(loc='best')

    return fig, axs

# MCMC
def plot_traces(list_of_tr, names=None):
    if names is None:
        names = np.arange(len(list_of_tr))
    # print (names)
    # Plot the variables
    # Initialize data frame for later analysis
    from scipy.stats.mstats import mquantiles
    # tr_len = mcmc.trace(stoch[0])[:, None].shape[0]
    # Height is determined by number of variables
    h = 3 * len(list_of_tr)
    fig, axs = plt.subplots(len(list_of_tr), 2, figsize=(16, h))

    for i, tr in enumerate(list_of_tr):
        try:
            # print str(tr)
            # Prepare Values
            tr_val = tr
            quants = mquantiles(tr_val, prob=[0.025, 0.25, 0.5, 0.75, 0.975])

            # Plot
            ### Left: Histogram
            h = axs[i, 0].hist(tr_val, histtype='stepfilled',
                               bins=50, label=str(names[i]), alpha=0.9)
            m = max(h[0])
            # m=1
            axs[i, 0].fill_betweenx([0, m + 5], quants[0], quants[-1], color='lightgreen')
            axs[i, 0].fill_betweenx([0, m + 5], quants[1], quants[-2], color='darkgreen')
            # axs[i, 0].set_title("{var} = {value:.4f} [{lci:.4f}, {hci:.4f}] ({tmp:.4f})".format(var=str(names[i]),
            #                                                                                tmp=tr_val.mean(),
            #                                                                                value=quants[2],
            #                                                                                lci=quants[0],
            #                                                                                hci=quants[-1]))
            axs[i, 0].set_title("{var} = {value:.4f} [{lci:.4f}, {hci:.4f}] ({tmp:.4f})".format(var=names[i],
                                                                                tmp=tr_val.mean(),
                                                                                value=quants[2],
                                                                                lci=quants[0],
                                                                                hci=quants[-1]))
            # h = axs[i, 0].hist(tr_val, histtype='stepfilled',
            #                    bins=50, label=str(names[i]), alpha=0.6)

            axs[i, 0].set_xticks(axs[i, 0].get_xticks()[::2])  # X 0
            axs[i, 0].set_yticks(axs[i, 0].get_yticks()[::2])  # Y 0
            axs[i, 0].set_ylim([0, m + 5])
            # print axs[i,1].get_yticks()[::2]
            ### Right: Trace
            l = range(len(tr_val))
            axs[i, 1].plot(l, tr_val)
            axs[i, 1].set_xticks(axs[i, 1].get_xticks()[::2])  # X 1
            axs[i, 1].set_yticks(axs[i, 1].get_yticks()[::2])  # Y 1
            axs[i, 1].set_xlim([0, l[-1]])

        except Exception as e:
            print (e)
    # fig.set_tight_layout
    fig.suptitle("Posterior Distribution <|> Convergence")
    return fig, axs

def plot_monthly(x_axis, data, model):
    fig, axs = plt.subplots(1, 2, figsize=(18, 8))
    # print(z.shape)
    labels = ['0-1', '1-21', '21+']
    ax = axs[0]
    for i, l in enumerate(labels):
        ax.plot(x_axis, model[i, :], label=l)
        ax.scatter(x_axis, data[i, :])
    ax.legend()

    ax = axs[1]
    ax.plot(x_axis, model.sum(axis=0), color='grey')
    ax.scatter(x_axis, data.sum(axis=0), c='grey')

    return fig, axs


text:
#-----------------------------------------------------------------------
# Copyright (C) 2020, All rights reserved
#
# Jonathan L. Hodges
#
#-----------------------------------------------------------------------
#=======================================================================
#
# DESCRIPTION:
# This software is part of a python library to assist in developing and
# analyzing simulation results from Fire Dynamics Simulator (FDS).
# FDS is an open source software package developed by NIST. The source
# code is available at: https://github.com/firemodels/fds
#
#=======================================================================
# # IMPORTS
#=======================================================================
import pyfdstools as fds
import os
from collections import defaultdict
import numpy as np
import matplotlib.pyplot as plt
import bpy

def exampleImportFile(fdsPath=None):
    if fdsPath is None:
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        fdsPath = os.path.join(systemPath, "examples", "case001.fds")
    fdsFile = fds.fdsFileOperations()
    print("\tStarting to import file from %s"%(fdsPath))
    fdsFile.importFile(fdsPath)
    print("\tFinished importing file.")
    return fdsFile

def exampleSaveFile(file=None, outdir=None, outname=None):
    if file is None:
        file = exampleImportFile()
    if outdir is None:
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        outdir = os.path.join(systemPath, "generated")
    if outname is None:
        location = os.path.join(outdir, '%s.fds'%(file.head['ID']['CHID']))
    print("\tStarting to save model.")
    file.saveModel(1, location, allowMeshSplitting=False)
    print("\tFinished saving model.")
    return location

def exampleErrorCalculation():
    data, keys = fds.readErrorTable()

    errorvalues = fds.calculatePercentile([100, 200, 300, 400, 500, 600], 'Plume Temperature', 0.95, fdsVersion='6.7.1')
    fig, ax1 = fds.plotPercentile(500, 'Plume Temperature', fdsVersion='6.7.1')
    return errorvalues

def exampleReadSlcf2dResults(resultDir=None, chid=None,
                             fdsQuantities = ['TEMPERATURE'],
                             tStart=0, tEnd=120):
    if (resultDir is None) and (chid is None):
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        chid = "case001"
        resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))

    datas = defaultdict(bool)

    for qty in fdsQuantities:
        grid, data, times = fds.readSLCF2Ddata(chid, resultDir, qty)
        tStartInd = np.argwhere(times >= tStart)[0][0]
        tEndInd = np.argwhere(times <= tEnd)[-1][0]
        data_tavg = np.nanmean(data[:, :, :, tStartInd:tEndInd], axis=3)
        datas[qty] = data_tavg.copy()
    datas['GRID'] = grid
    datas['TIMES'] = times

    return datas

def exampleExtract2dFromSlcf2d(datas,
                               fdsQuantities = ['TEMPERATURE'],
                               fdsUnits = ['C'],
                               axis=1, value=2.5,
                               qnty_mn=20, qnty_mx=150):
    datas2D = defaultdict(bool)
    for qty, unit in zip(fdsQuantities, fdsUnits):
        x, z, data_slc = fds.findSliceLocation(datas['GRID'], datas[qty], axis, value)
        datas2D[qty] = data_slc
        fig = fds.plotSlice(x, z, data_slc, axis,
                            qnty_mn=qnty_mn, qnty_mx=qnty_mx,
                            clabel="%s (%s)"%(qty, unit),
                            cbarticks=[20, 40, 60, 80, 100, 120, 140])
    datas2D['X'] = x
    datas2D['Z'] = z

    return datas2D, fig

def exampleReadSlcf3dResults(resultDir=None, chid=None,
                             fdsQuantities = ['TEMPERATURE'],
                             tStart=0, tEnd=120):
    if (resultDir is None) and (chid is None):
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        chid = "case001"
        resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))

    datas = defaultdict(bool)

    for qty in fdsQuantities:
        grid, data, times = fds.readSLCF3Ddata(chid, resultDir, qty)
        tStartInd = np.argwhere(times >= tStart)[0][0]
        tEndInd = np.argwhere(times <= tEnd)[-1][0]
        data_tavg = np.nanmean(data[:, :, :, tStartInd:tEndInd], axis=3)
        datas[qty] = data_tavg.copy()
    datas['GRID'] = grid
    datas['TIMES'] = times

    return datas

def exampleExtract2dFromSlcf3d(datas,
                               fdsQuantities = ['TEMPERATURE'],
                               fdsUnits = ['C'],
                               axis=2, value=4.4,
                               qnty_mn=20, qnty_mx=150):
    datas2D = defaultdict(bool)
    for qty, unit in zip(fdsQuantities, fdsUnits):
        x, z, data_slc = fds.findSliceLocation(datas['GRID'], datas[qty], axis, value)
        datas2D[qty] = data_slc
        fig = fds.plotSlice(x, z, data_slc, axis,
                            qnty_mn=qnty_mn, qnty_mx=qnty_mx,
                            clabel="%s (%s)"%(qty, unit),
                            cbarticks=[20, 40, 60, 80, 100, 120, 140])
    datas2D['X'] = x
    datas2D['Z'] = z

    return datas2D, fig

def exampleImportBndf(resultDir=None, chid=None, fdsFile=None,
                       fdsQuantities = ['WALL TEMPERATURE'],
                       fdsUnits = ['C'],
                       tStart=0, tEnd=120,
                       axis=-2, value=4.4,
                       qnty_mn=20, qnty_mx=100):
    if (resultDir is None) and (chid is None):
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        chid = "case001"
        resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))

    if fdsFile is None: fdsFilePath = fds.getFileList(resultDir, chid, 'fds')[0]

    datas, times = fds.queryBndf(resultDir, chid, fdsFilePath, fdsQuantities, fdsUnits, axis, value)
    tStartInd = np.argwhere(times >= tStart)[0][0]
    tEndInd = np.argwhere(times <= tEnd)[-1][0]

    for qty, unit in zip(fdsQuantities, fdsUnits):
        for mesh in list(datas[qty].keys()):
            data = datas[qty]['DATA']
            x = datas[qty]['X']
            z = datas[qty]['Z']
        meanData = np.mean(data[:, :, tStartInd:tEndInd], axis=2)
        fig = fds.plotSlice(x, z, meanData, axis,
                            qnty_mn=qnty_mn, qnty_mx=qnty_mx,
                            clabel="%s (%s)"%(qty, unit),
                            cbarticks=[20, 40, 60, 80, 100, 120, 140])
    return datas, fig

def exampleExtractBndfMax(resultDir=None, chid=None, outDir=None,
                          fdsFile=None, smvFile=None, outputNamespace=None,
                          fdsQuantities = ['WALL TEMPERATURE'], fdsUnits = ['C'],
                          tStart=0, tEnd=120, tInt=1, tBand=3, orientations=[0],
                          axis=-2, value=4.4,
                          qnty_mn=20, qnty_mx=100,
                          yticks=[20, 50, 100, 150, 200, 250, 300, 350],
                          plot=True):
    if (resultDir is None) and (chid is None) and (outDir is None):
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        chid = "case001"
        resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))
        outDir = os.path.join(systemPath, "generated")

    if fdsFile is None: fdsFilePath = fds.getFileList(resultDir, chid, 'fds')[0]
    if smvFile is None: smvFilePath = fds.getFileList(resultDir, chid, 'smv')[0]
    if outputNamespace is None: outputNamespace = "%s_max_"%(chid)

    outputName = os.path.join(outDir, outputNamespace)

    datas = fds.extractMaxBndfValues(fdsFilePath, smvFilePath, resultDir, chid, fdsQuantities,
                    tStart=tStart, tEnd=tEnd, tInt=tInt, tBand=tBand, orientations=orientations)

    figs = []
    for qty in fdsQuantities:
        times = datas[qty]['TIMES']
        mPts = datas[qty]['DATA']
        names = datas[qty]['NAMES']
        outName = "%s%s"%(outputName, qty)
        fds.maxValueCSV(times, mPts, names, outName)
        if plot:
            fig = fds.maxValuePlot(times, mPts, names, outName, vName=qty, yticks=yticks)
            figs.append(fig)
    return datas, figs

def example2dSliceToCsv(resultDir=None, outDir=None, chid=None,
                        quantity=None, unit=None, axis=None, value=None,
                        time=None, dt=None):
    if (resultDir is None) and (chid is None) and (outDir is None):
        try:
            import bpy
            systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
        except:
            systemPath = os.path.dirname(os.path.abspath(__file__))
        chid = "case001"
        resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))
        outDir = os.path.join(systemPath, "generated")
    try:
        os.mkdir(outDir)
    except:
        pass
    data = fds.query2dAxisValue(resultDir, chid, quantity, axis, value, time=time, dt=dt)
    fds.renderSliceCsvs(data, chid, outDir)
    fig, ax = fds.plotSlice(data['x'], data['z'], data['datas'][:, :, -1], axis,
                        clabel="%s (%s)"%(quantity, unit))
    fig.savefig(os.path.join(outDir, '%s_%s_%0.0f_%0.4f_final_frame.png'%(chid, quantity, axis, value)))
    return data, fig

def runExamples():
    try:
        import bpy
        systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
    except:
        systemPath = os.path.dirname(os.path.abspath(__file__))
    exampleInputFdsFile = os.path.join(systemPath, "examples", "case001.fds")
    exampleOutputDir = os.path.join(systemPath, "generated")
    chid = "case001"
    resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))

    print("Importing model example", flush=True)
    file = exampleImportFile(exampleInputFdsFile)

    print("Saving model example", flush=True)
    exampleSaveFile(file=file, outdir=exampleOutputDir)

    print("Plotting error example", flush=True)
    exampleErrorCalculation()

    print("Importing SLCF2D results example", flush=True)
    datas = exampleReadSlcf2dResults(resultDir=resultDir, chid=chid)

    print("Extracting slice from SLCF2D results example", flush=True)
    datas2D, fig = exampleExtract2dFromSlcf2d(datas)

    print("Importing SLCF3D results example", flush=True)
    datas = exampleReadSlcf3dResults(resultDir=resultDir, chid=chid)

    print("Extracting 2-D slice from SLCF3D results example", flush=True)
    datas2D, fig = exampleExtract2dFromSlcf3d(datas)

    print("Importing BNDF results example", flush=True)
    datas, fig = exampleImportBndf(resultDir=resultDir, chid=chid)

    print("Extracting max value from BNDF results example", flush=True)
    datas, figs = exampleExtractBndfMax(resultDir=resultDir, chid=chid, outDir=exampleOutputDir)

    print("Rendering 2d slice to csv example.", flush=True)
    datas = example2dSliceToCsv(resultDir=resultDir, chid=chid, outDir=exampleOutputDir,
                                axis=1, value=2.45, time=30, dt=60, quantity='TEMPERATURE', unit='C')

    return systemPath, exampleInputFdsFile, exampleOutputDir, chid, resultDir

def create_custom_mesh(objname, px, py, pz):

    # Define arrays for holding data
    myvertex = []
    myfaces = []

    # Create all Vertices

    # vertex 0
    mypoint = [(-1.0, -1.0, 0.0)]
    myvertex.extend(mypoint)

    # vertex 1
    mypoint = [(1.0, -1.0, 0.0)]
    myvertex.extend(mypoint)

    # vertex 2
    mypoint = [(-1.0, 1.0, 0.0)]
    myvertex.extend(mypoint)

    # vertex 3
    mypoint = [(1.0, 1.0, 0.0)]
    myvertex.extend(mypoint)

    # -------------------------------------
    # Create all Faces
    # -------------------------------------
    myface = [(0, 1, 3, 2)]
    myfaces.extend(myface)

    mymesh = bpy.data.meshes.new(objname)

    myobject = bpy.data.objects.new(objname, mymesh)

    bpy.context.scene.collection.objects.link(myobject)

    # Generate mesh data
    mymesh.from_pydata(myvertex, [], myfaces)
    # Calculate the edges
    mymesh.update(calc_edges=True)

    # Set Location
    myobject.location.x = px
    myobject.location.y = py
    myobject.location.z = pz

    return myobject

def create_custom_mesh2(objname, x0, x1, y0, y1, z0, z1, data):

    # Define arrays for holding data
    myvertex = []
    myfaces = []

    # Create all Vertices

    # vertex 0
    mypoint = [(x0, y0, z0)]
    myvertex.extend(mypoint)

    # vertex 1
    mypoint = [(x1, y0, z0)]
    myvertex.extend(mypoint)

    # vertex 2
    mypoint = [(x0, y1, z1)]
    myvertex.extend(mypoint)

    # vertex 3
    mypoint = [(x1, y1, z1)]
    myvertex.extend(mypoint)

    # -------------------------------------
    # Create all Faces
    # -------------------------------------
    myface = [(0, 1, 3, 2)]
    myfaces.extend(myface)

    print("A", objname)
    mymesh = bpy.data.meshes.new(objname)
    print("B")
    myobject = bpy.data.objects.new(objname, mymesh)
    print("C")
    bpy.context.scene.collection.objects.link(myobject)
    print("D")
    # Generate mesh data
    mymesh.from_pydata(myvertex, [], myfaces)
    # Calculate the edges
    print("E")
    mymesh.update(calc_edges=True)

    # Set Location
    #myobject.location.x = px
    #myobject.location.y = py
    #myobject.location.z = pz
    print("F")
    myobject.location.x = 0
    myobject.location.y = 0
    myobject.location.z = 0

    tex = bpy.data.textures.new('%s-TEXTURE'%(objname), 'IMAGE')
    print(data.shape)
    mn = np.nanmin(data)
    mx = np.nanmax(data)
    d = (data - mn)/(mx - mn)
    d2 = np.zeros((d.shape[0], d.shape[1], 4))
    d2[:, :, 0] = d
    d2[:, :, 1] = d
    d2[:, :, 2] = d
    d2[:, :, 3] = 1
    print(d2.shape, d2.flatten().shape)

    image = bpy.data.images.new("%s-IMAGE"%(objname), width=d.shape[0], height=d.shape[1])
    image.pixels = d2.flatten()
    tex.image = image

    mat = bpy.data.materials.new('%s-DATA'%(objname))
    mat.texture_paint_images = [image]
    #mat.texture_slots.add()
    #ts = mat.texture_slots[0]
    #ts.texture = tex
    #ts.texture_coords = 'UV'
    #ts.uv_layer = 'default'
    #mtex.texture_coords = 'UV'

    myobject.data.materials = mat

    return myobject

if __name__ == '__main__':
    QTPLUGINPATH = os.getenv('QT_PLUGIN_PATH')
    print(QTPLUGINPATH)

    try:
        import bpy
        systemPath = os.path.dirname(bpy.context.space_data.text.filepath)
    except:
        systemPath = os.path.dirname(os.path.abspath(__file__))
    exampleInputFdsFile = os.path.join(systemPath, "examples", "case001.fds")
    exampleOutputDir = os.path.join(systemPath, "generated")
    chid = "case001"
    resultDir = os.path.join(systemPath, "examples", "%s.zip"%(chid))

    #datas, fig = exampleImportBndf(resultDir=resultDir, chid=chid)
    #print(datas['WALL TEMPERATURE'].keys(), flush=True)
    #print(datas.keys(), flush=True)
    datas, figs = exampleExtractBndfMax(resultDir=resultDir, chid=chid, outDir=exampleOutputDir, plot=False)
    times = datas['WALL TEMPERATURE']['TIMES']
    names = datas['WALL TEMPERATURE']['NAMES']
    ds = datas['WALL TEMPERATURE']['DATA']

    import pandas as pd
    print(ds.shape, flush=True)
    print(len(times), flush=True)
    print(len(names), flush=True)
    tmp = pd.DataFrame(ds, index=times, columns=names)
    print(tmp, flush=True)

    qty = 'WALL TEMPERATURE'
    for key in names:
        bpy.context.scene.collection.objects[key]['%s TIME'%(qty)] = tmp.index
        bpy.context.scene.collection.objects[key]['%s MAX'%(qty)] = tmp[key]

    systemPath, exampleInputFdsFile, exampleOutputDir, chid, resultDir = runExamples()

    print(bpy.context.scene.collection.objects['COUCH-BACK-XPOS']['WALL TEMPERATURE MAX'][:])
    print(bpy.context.scene.collection.objects['COUCH-BACK-XPOS']['WALL TEMPERATURE TIME'][:])

    '''
    qty = 'TEMPERATURE'
    axis = 2
    value = 4.4
    obst = 'Centerline'

    datas = exampleReadSlcf3dResults(resultDir=resultDir, chid=chid)
    datas2D, fig = exampleExtract2dFromSlcf3d(datas, axis=axis, value=value)

    TEMPERATURE = datas2D['TEMPERATURE']
    X = datas2D['X']
    Z = datas2D['Z']

    x0 = X[0, 0]
    x1 = X[-1, -1]
    y0 = value
    y1 = value
    z0 = Z[0, 0]
    z1 = Z[-1, -1]
    corners = [(x0, y0, z0), (x1, y0, z0), (x1, y0, z1), (x0, y0, z1)]
    #create_custom_mesh("Awesome_object", curloc[0], curloc[1], 0)
    objname = "%s_%s"%(obst, qty)
    print(objname)
    print(x0, x1, y0, y1, z0, z1)
    #myobj = create_custom_mesh2(objname, x0, x1, y0, y1, z0, z1, TEMPERATURE)
    #print(myobj.__dict__)
    #x, z, data_slc = fds.findSliceLocation(datas['GRID'], datas[qty], axis, value)
    '''

    '''
    print(x[:, 0], x[0, :])
    #corners = [(x[0, 0], value, z[0, 0])]
    #qty = 'TEMPERATURE'
    #grid, data, times = fds.readSLCF3Ddata(chid, resultDir, qty)

    #corners =

    print(grid.shape, data.shape, times.shape)

    d = np.nanmean(datas['TEMPERATURE'], axis=2)
    import bpy
    context = bpy.context

    curloc = bpy.context.scene.cursor.location

    create_custom_mesh("Awesome_object", curloc[0], curloc[1], 0)
    '''

max_stars_repo_name: johodges/pyfdstools

max_issues_repo_name: johodges/pyfdstools

max_forks_repo_name: johodges/pyfdstools


text:
from __future__ import division

import numpy as np
import pylab as plt
from scipy import stats

def hist(ax, stars_characteristics, predictions):
	"""

	"""

	binary_stars = stars_characteristics[:,0]
	idbin = np.where(binary_stars == 1)

	all_stars = stars_characteristics[idbin, 1].flatten() / .12 # Convert from Euclid small px to mas
	_, binst, _ = ax.hist(all_stars, 100, histtype='stepfilled', color="r", alpha=0.8)

	observed = np.where(np.logical_and(binary_stars == 1, predictions))
	obs_stars = stars_characteristics[observed, 1].flatten() / .12 # Convert from Euclid small px to mas
	ax.hist(obs_stars, binst, histtype='stepfilled', color="darkgrey", alpha=0.8)

	binst = stats.mstats.mquantiles(all_stars, [0.,0.2,0.4,0.5,0.6,0.7,0.8,0.9,0.99,1.])
	#binst = 6
	nt, binst = np.histogram(all_stars, binst)
	nd, _ = np.histogram(obs_stars, bins=binst)
	xbin = [(binst[i] + binst[i+1])/2. for i in range(len(binst) - 1)]
	ybin = nd/nt
	#xbin = np.concatenate([[all_stars.min()], xbin])
	#ybin = np.concatenate([[ybin[0]], ybin])

	ax2 = ax.twinx()
	ax2.set_ylabel(r'$\mathrm{Completeness}$')
	ax2.plot(xbin, ybin, c='k', lw=2)
	ax2.set_ylim([0, 1])

	plt.setp(ax.get_yticklabels(), visible=False)

	ax2.grid(True)
	ax.set_xlabel(r"$\mathrm{Angular\ separation\ [mas]}$")
	ax.set_ylabel(r"$\mathrm{Occurence\ of\ binaries\ (arbit.\ units)}$")

	ax.set_xlim([all_stars.min(), (all_stars.max())])

max_stars_repo_name: kuntzer/binfind

max_issues_repo_name: kuntzer/binfind

max_forks_repo_name: kuntzer/binfind