text:
import numpy as np
from environments.gridworld import GridWorld
from utils.util import softmax

GRID_HEIGHT = 4
GRID_WIDTH = 4
TERMINAL_STATES = [(0, 0), (GRID_HEIGHT-1, GRID_WIDTH-1)]
DISCOUNT_RATE = 1.0
THETA_1 = 0.0001
THETA_2 = 0.0001

MAX_EPISODES = 5000

# 정책 반복 클래스
class PolicyIteration:
    def __init__(self, env):
        self.env = env

        self.max_iteration = MAX_EPISODES

        self.terminal_states = [(0, 0), (4, 4)]

        self.state_values = None
        self.policy = self.generate_random_policy()

        # 정책 평가
        self.delta = 0.0

        # 정책 평가의 역
        self.theta = 0.001

    # 모든 상태에서 수행 가능한 행동에 맞춰 임의의 정책을 생성함
    # 초기에 각 행동의 선택 확률은 모두 같음
    def generate_random_policy(self):
        policy = dict()

        for i in range(GRID_HEIGHT):
            for j in range(GRID_WIDTH):
                actions = []
                action_probs = []
                for action in range(self.env.action_space.num_actions):
                    actions.append(action)
                    action_probs.append(0.25)
                policy[(i, j)] = (actions, action_probs)

        return policy

    # 정책 평가 함수
    def policy_evaluation(self):
        # 상태-가치 함수 초기화
        state_values = np.zeros((GRID_HEIGHT, GRID_WIDTH))

        # 가치 함수의 값들이 수렴할 때까지 반복
        iter_num = 0
        while True:
            old_state_values = state_values.copy()

            for i in range(GRID_HEIGHT):
                for j in range(GRID_WIDTH):
                    if (i, j) in TERMINAL_STATES:
                        state_values[i][j] = 0.0
                    else:
                        values = []
                        for action in self.env.action_space.ACTIONS:
                            (next_i, next_j), reward, prob = env.get_state_action_probability(state=(i, j), action=action)

                            _, action_probs = self.policy[(i, j)]

                            # Bellman-Equation, 벨만 방정식 적용
                            values.append(
                                action_probs[action] * prob * (reward + DISCOUNT_RATE * old_state_values[next_i, next_j])
                            )

                        state_values[i][j] = np.sum(values)

            iter_num += 1

            # 갱신되는 값이 THETA_1(=0.0001)을 기준으로 수렴하는지 판정
            max_delta_value = abs(old_state_values - state_values).max()
            if max_delta_value < THETA_1:
                break

        self.state_values = state_values

        return iter_num

    # 정책 개선 함수
    def policy_improvement(self):
        new_policy = dict()

        is_policy_stable = True

        # 행동-가치 함수 생성
        for i in range(GRID_HEIGHT):
            for j in range(GRID_WIDTH):
                if (i, j) in TERMINAL_STATES:
                    actions = []
                    action_probs = []
                    for action in range(self.env.action_space.num_actions):
                        actions.append(action)
                        action_probs.append(0.25)
                    new_policy[(i, j)] = (actions, action_probs)
                else:
                    actions = []
                    q_func = []
                    for action in self.env.action_space.ACTIONS:
                        actions.append(action)
                        (next_i, next_j), reward, prob = env.get_state_action_probability(state=(i, j), action=action)
                        q_func.append(
                            prob * (reward + DISCOUNT_RATE * self.state_values[next_i, next_j])
                        )

                    new_policy[(i, j)] = (actions, softmax(q_func))

        error = 0.0
        for i in range(GRID_HEIGHT):
            for j in range(GRID_WIDTH):
                error += np.sum(np.absolute(np.array(self.policy[(i, j)][1]) - np.array(new_policy[(i, j)][1])))

        if error > THETA_2:
            is_policy_stable = False

        self.policy = new_policy

        return is_policy_stable, error

    # 정책 반복 함수
    def start_iteration(self):
        iter_num = 0

        # 정책의 안정성 검증
        is_policy_stable = False

        print("[[[ 정책 반복 시작! ]]]")
        while not is_policy_stable and iter_num < self.max_iteration:
            iter_num_policy_evaluation = self.policy_evaluation()
            print("*** 정책 평가 [수렴까지 누적 반복 횟수: {0}] ***".format(iter_num_policy_evaluation))

            is_policy_stable, error = self.policy_improvement()
            print("*** 정책 개선 [에러 값: {0:7.5f}] ***".format(error))

            iter_num += 1
            print("*** 정책의 안정(Stable) 여부: {0}, 반복 횟수: {1} ***\n".format(
                is_policy_stable,
                iter_num
            ))
        print("[[[ 정책 반복 종료! ]]]\n\n")

        return self.state_values, self.policy

if __name__ == '__main__':
    # 그리드 월드 환경 객체 생성
    env = GridWorld(
        height=GRID_HEIGHT,
        width=GRID_WIDTH,
        start_state=(0, 0),
        terminal_states=TERMINAL_STATES,
        transition_reward=-1.0,
        terminal_reward=-1.0,
        outward_reward=-1.0
    )
    env.reset()

    PI = PolicyIteration(env)
    PI.start_iteration()
    print(PI.state_values, end="\n\n")

    for i in range(GRID_HEIGHT):
        for j in range(GRID_WIDTH):
            print(
                i, j,
                ": UP, DOWN, LEFT, RIGHT",
                PI.policy[(i, j)][1]
            )
        print()


text:
#
# Running the validation score
# for the steering models
#
# (c) Yongyang Nie
# 2018, All Rights Reserved
#

from i3d import Inception3D
import configs
from os import path
import pandas as pd
import numpy as np
import helper
import math
import time
import communication

def validation_score(model_path, write_output=False):

    model = Inception3D(input_shape=(configs.LENGTH, configs.IMG_HEIGHT, configs.IMG_WIDTH, configs.CHANNELS),
                        weights_path=model_path)

    # steerings and images
    steering_labels = path.join(configs.VAL_DIR, 'labels.csv')

    df_truth = pd.read_csv(steering_labels, usecols=['frame_id', 'steering_angle'], index_col=None)

    esum = 0
    inputs = []
    predictions = []
    start_time = time.time()

    for i in range(configs.LENGTH):
        file = configs.VAL_DIR + "center/" + str(df_truth['frame_id'].loc[i]) + ".jpg"
        img = helper.load_image(file)
        inputs.append(img)

    # Run through all images
    for i in range(configs.LENGTH, len(df_truth)):

        img = helper.load_image(configs.VAL_DIR + "center/" + str(df_truth['frame_id'].loc[i]) + ".jpg")
        inputs.pop(0)
        inputs.append(img)
        prediction = model.model.predict(np.array([inputs]))[0]
        prediction = prediction[0]
        actual_steers = df_truth['steering_angle'].loc[i]
        e = (actual_steers - prediction) ** 2
        esum += e

        predictions.append(prediction)

        if len(predictions) % 1000 == 0:
            print('.')

    print("time per step: %s seconds" % ((time.time() - start_time) / len(predictions)))

    if write_output:
        print("Writing predictions...")
        pd.DataFrame({"steering_angle": predictions}).to_csv('./result.csv', index=False, header=True)
        print("Done!")
    else:
        print("Not writing outputs")
        print("Done")

    return math.sqrt(esum / len(predictions))

if __name__ == "__main__":

    print("Validating...")
    score = validation_score('i3d_rgb_64_v10.h5')
    print("score: " + str(score))

    communication.notify_validation_completion(score, 'i3d_rgb_64_v10.h5')

ext: py

max_stars_repo_name: sigmaai/behavioral-cloning

max_issues_repo_name: sigmaai/behavioral-cloning

max_forks_repo_name: sigmaai/behavioral-cloning

reason: import numpy


text:
'''
Preprocess pix3d data.
author: ynie
date: Sep, 2019
'''

import sys
sys.path.append('.')
from configs.pix3d_config import Config
import json
from multiprocessing import Pool
import pickle
from libs.tools import read_obj, sample_pnts_from_obj, normalize_to_unit_square
import os
from PIL import Image
import numpy as np

def parse_sample(sample_info):
    sample_id = sample_info[0]
    print('Processing sample: %d.' % (sample_id))
    sample = sample_info[1]
    is_flip = sample_info[2]

    class_id = config.classnames.index(sample['category'])
    save_path = os.path.join(config.train_test_data_path, str(sample_id) + '.pkl')

    if os.path.exists(save_path):
        return None

    data = {}
    data['sample_id'] = sample_id

    '''get single object image'''
    img = np.array(Image.open(os.path.join(config.metadata_path, sample['img'])).convert('RGB'))
    img = img[sample['bbox'][1]:sample['bbox'][3], sample['bbox'][0]:sample['bbox'][2]]
    data['img'] = img

    '''get object model'''
    model_path = os.path.join(config.metadata_path, sample['model'])
    obj_data = read_obj(model_path, ['v', 'f'])
    sampled_points = sample_pnts_from_obj(obj_data, 10000, mode='random')
    sampled_points = normalize_to_unit_square(sampled_points)[0]

    data['gt_3dpoints'] = sampled_points

    '''get object category'''
    data['class_id'] = class_id

    pickle.dump(data, save_path)

    if is_flip:
        save_path = os.path.join(config.train_test_data_path, str(sample_id) + '_flipped.pkl')
        if os.path.exists(save_path):
            return None

        data_flip = {}
        data_flip['sample_id'] = sample_id
        img_flip = np.array(Image.open(os.path.join(config.metadata_path, sample['img'])).convert('RGB').transpose(Image.FLIP_LEFT_RIGHT))
        x1 = sample['img_size'][0] - 1 - sample['bbox'][2]
        x2 = sample['img_size'][0] - 1 - sample['bbox'][0]
        y1 = sample['bbox'][1]
        y2 = sample['bbox'][3]

        img_flip = img_flip[y1:y2, x1:x2]
        data_flip['img'] = img_flip
        data_flip['gt_3dpoints'] = sampled_points
        data_flip['class_id'] = class_id
        pickle.dump(data_flip, save_path)

if __name__ == '__main__':
    config = Config('pix3d')
    with open(config.metadata_file, 'r') as file:
        metadata = json.load(file)

    with open(config.train_split, 'r') as file:
        train_split = json.load(file)

    with open(config.test_split, 'r') as file:
        test_split = json.load(file)

    train_ids = [int(os.path.basename(file).split('.')[0]) for file in train_split if 'flipped' not in file]
    test_ids = [int(os.path.basename(file).split('.')[0]) for file in test_split if 'flipped' not in file]

    train_sample_info = [(sample_id, metadata[sample_id], True) for sample_id in train_ids]
    test_sample_info = [(sample_id, metadata[sample_id], False) for sample_id in test_ids]

    p = Pool(processes=1)
    p.map(parse_sample, train_sample_info)
    p.close()
    p.join()

    p = Pool(processes=1)
    p.map(parse_sample, test_sample_info)
    p.close()
    p.join()


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Tue Oct  1 17:05:20 2019

@author: jorge.guerra
"""
import numpy as np
import scipy.linalg as dsl
import matplotlib.pyplot as plt
import scipy.optimize as opt
import scipy.sparse.linalg as spl
import computeEulerEquationsLogPLogT as tendency
from computeResidualViscCoeffs import computeResidualViscCoeffs

def pause():
    input("Press the <ENTER> key to continue...")

def computeIterativeSolveNL(PHYS, REFS, REFG, DX, DZ, SOLT, INIT, udex, wdex, pdex, tdex, botdex, topdex, wbdex, sysDex, DynSGS):
       dHdX = REFS[6]
       lastSol = SOLT[:,0]

       def computeRHSUpdate(sol):
              # Update the BC
              sol[wbdex] = dHdX * (INIT[botdex] + sol[botdex])
              # Function evaluation
              fields, U, RdT = tendency.computePrepareFields(PHYS, REFS, sol, INIT, udex, wdex, pdex, tdex)
              rhs = tendency.computeEulerEquationsLogPLogT_NL(PHYS, REFS, REFG, fields, U, RdT)
              rhs += tendency.computeRayleighTendency(REFG, fields)

              if DynSGS:
                     RESCF = computeResidualViscCoeffs(sol, rhs, DX, DZ, udex, wdex, pdex, tdex)
                     rhs += tendency.computeDynSGSTendency(RESCF, REFS, fields, udex, wdex, pdex, tdex, botdex, topdex)

              return rhs

       def computeJacVecUpdate(sol, vec):
              # Prepare the Jacobian
              fields, U, RdT = tendency.computePrepareFields(PHYS, REFS, sol, INIT, udex, wdex, pdex, tdex)
              DOPS = tendency.computeJacobianMatrixLogPLogT(PHYS, REFS, REFG, fields, U, RdT, botdex, topdex)
              # Compute the product
              jv = tendency.computeJacobianVectorProduct(DOPS, REFG, vec, udex, wdex, pdex, tdex)

              return jv

       # Solve for nonlinear equilibrium (default krylov)
       #'''
       jac_options = {'method':'gmres','inner_maxiter':1000,'outer_k':3}
       sol, info = opt.nonlin.nonlin_solve(computeRHSUpdate, lastSol,
                                  jacobian=opt.nonlin.KrylovJacobian(**jac_options),
                                  iter=5, verbose=True,
                                  maxiter=100,
                                  line_search='armijo',
                                  full_output=True,
                                  raise_exception=False)

       # Solve for nonlinear equilibrium (modified direct krylov)
       '''
       jac_options = {'jdv':computeJacVecUpdate, \
                      'method':'gmres', \
                      'inner_maxiter':2000, \
                      'outer_k':2}
       sol, info = opt.nonlin.nonlin_solve(computeRHSUpdate, lastSol,
                                  jacobian=KrylovDirectJacobian(**jac_options),
                                  iter=5, verbose=True,
                                  maxiter=100,
                                  line_search='armijo',
                                  full_output=True,
                                  raise_exception=False)
       '''
       print('NL solver exit on: ', info)

       return sol

#------------------------------------------------------------------------------
# Iterative/Krylov un-approximated Jacobians
#------------------------------------------------------------------------------

class KrylovDirectJacobian(opt.nonlin.Jacobian):
    r"""
    Find a root of a function, using Krylov approximation for inverse Jacobian.

    This method is suitable for solving large-scale problems.

    Parameters
    ----------
    %(params_basic)s
    rdiff : float, optional
        Relative step size to use in numerical differentiation.
    method : {'lgmres', 'gmres', 'bicgstab', 'cgs', 'minres'} or function
        Krylov method to use to approximate the Jacobian.
        Can be a string, or a function implementing the same interface as
        the iterative solvers in `scipy.sparse.linalg`.

        The default is `scipy.sparse.linalg.lgmres`.
    inner_M : LinearOperator or InverseJacobian
        Preconditioner for the inner Krylov iteration.
        Note that you can use also inverse Jacobians as (adaptive)
        preconditioners. For example,

        >>> from scipy.optimize.nonlin import BroydenFirst, KrylovJacobian
        >>> from scipy.optimize.nonlin import InverseJacobian
        >>> jac = BroydenFirst()
        >>> kjac = KrylovJacobian(inner_M=InverseJacobian(jac))

        If the preconditioner has a method named 'update', it will be called
        as ``update(x, f)`` after each nonlinear step, with ``x`` giving
        the current point, and ``f`` the current function value.
    inner_tol, inner_maxiter, ...
        Parameters to pass on to the \"inner\" Krylov solver.
        See `scipy.sparse.linalg.gmres` for details.
    outer_k : int, optional
        Size of the subspace kept across LGMRES nonlinear iterations.
        See `scipy.sparse.linalg.lgmres` for details.
    %(params_extra)s

    See Also
    --------
    scipy.sparse.linalg.gmres
    scipy.sparse.linalg.lgmres

    Notes
    -----
    This function implements a Newton-Krylov solver. The basic idea is
    to compute the inverse of the Jacobian with an iterative Krylov
    method. These methods require only evaluating the Jacobian-vector
    products, which are NOT approximated by evaluated directly and
    passed to the "matvec" member function.

    Due to the use of iterative matrix inverses, these methods can
    deal with large nonlinear problems.

    SciPy's `scipy.sparse.linalg` module offers a selection of Krylov
    solvers to choose from. The default here is `lgmres`, which is a
    variant of restarted GMRES iteration that reuses some of the
    information obtained in the previous Newton steps to invert
    Jacobians in subsequent steps.

    For a review on Newton-Krylov methods, see for example [1]_,
    and for the LGMRES sparse inverse method, see [2]_.

    References
    ----------
    .. [1] D.A. Knoll and D.E. Keyes, J. Comp. Phys. 193, 357 (2004).
           :doi:`10.1016/j.jcp.2003.08.010`
    .. [2] A.H. Baker and E.R. Jessup and T. Manteuffel,
           SIAM J. Matrix Anal. Appl. 26, 962 (2005).
           :doi:`10.1137/S0895479803422014`

    """

    def __init__(self, rdiff=None, jdv=None, method='lgmres', inner_maxiter=20,
                 inner_M=None, outer_k=10, **kw):
        self.rdiff = rdiff
        self.jac_vec = jdv
        self.preconditioner = inner_M
        self.method = dict(
            bicgstab = spl.bicgstab,
            gmres = spl.gmres,
            lgmres = spl.lgmres,
            cgs = spl.cgs,
            minres = spl.minres,
            ).get(method, method)

        self.method_kw = dict(maxiter=inner_maxiter, M=self.preconditioner)

        if self.method is spl.gmres:
            # Replace GMRES's outer iteration with Newton steps
            self.method_kw['restrt'] = inner_maxiter
            self.method_kw['maxiter'] = 1
            self.method_kw.setdefault('atol', 0)
        elif self.method is spl.gcrotmk:
            self.method_kw.setdefault('atol', 0)
        elif self.method is spl.lgmres:
            self.method_kw['outer_k'] = outer_k
            # Replace LGMRES's outer iteration with Newton steps
            self.method_kw['maxiter'] = 1
            # Carry LGMRES's `outer_v` vectors across nonlinear iterations
            self.method_kw.setdefault('outer_v', [])
            self.method_kw.setdefault('prepend_outer_v', True)
            # But don't carry the corresponding Jacobian*v products, in case
            # the Jacobian changes a lot in the nonlinear step
            #
            # XXX: some trust-region inspired ideas might be more efficient...
            #      See eg. Brown & Saad. But needs to be implemented separately
            #      since it's not an inexact Newton method.
            self.method_kw.setdefault('store_outer_Av', False)
            self.method_kw.setdefault('atol', 0)

        for key, value in kw.items():
            if not key.startswith('inner_'):
                raise ValueError("Unknown parameter %s" % key)
            self.method_kw[key[6:]] = value

    def _update_diff_step(self):
        mx = abs(self.x0).max()
        mf = abs(self.f0).max()
        self.omega = self.rdiff * max(1, mx) / max(1, mf)

    def matvec(self, v):
        nv = dsl.norm(v)
        if nv == 0:
            return 0*v
        # Compute the updated Jacobian-vector product
        #r1 = self.jac_vec(self.x0, v)
        sc = 0.5*self.omega / nv
        #r = (self.func(self.x0 + sc*v) - self.f0) / sc
        r2 = (self.func(self.x0 + sc*v) - self.func(self.x0 - sc*v)) / sc
        '''
        plt.figure(figsize=(12.0, 8.0))
        varL = (int) (len(v) / 4)
        udex = range(0, varL)
        wdex = range(varL, 2*varL)
        pdex = range(2*varL, 3*varL)
        tdex = range(3*varL, len(v))
        plt.plot(r1[tdex], 'k-'); #plt.xscale('log'); plt.yscale('linear'); #plt.ylim(1.0E-6, 1.0E-1)
        plt.plot(r2[tdex], 'b:'); #plt.xscale('log'); plt.yscale('linear'); #plt.ylim(1.0E-6, 1.0E-1)
        plt.show()
        pause()
        '''
        r = r2

        if not np.all(np.isfinite(r)) and np.all(np.isfinite(v)):
            raise ValueError('Function returned non-finite results')
        return r

    def solve(self, rhs, tol=0):
        if 'tol' in self.method_kw:
            sol, info = self.method(self.op, rhs, **self.method_kw)
        else:
            sol, info = self.method(self.op, rhs, tol=tol, **self.method_kw)
        return sol

    def update(self, x, f):
        self.x0 = x
        self.f0 = f
        self._update_diff_step()

    def setup(self, x, f, func):
        opt.nonlin.Jacobian.setup(self, x, f, func)
        self.x0 = x
        self.f0 = f
        self.op = spl.aslinearoperator(self)

        if self.rdiff is None:
            self.rdiff = np.finfo(x.dtype).eps ** (1./2)

        self._update_diff_step()

max_stars_repo_name: jeguerra/nonlinearMtnWavesSolver

max_issues_repo_name: jeguerra/nonlinearMtnWavesSolver

max_forks_repo_name: jeguerra/nonlinearMtnWavesSolver

reason: import numpy,import scipy,from scipy


text:
#!/usr/bin/env python

import pickle
import InstrumentDriver
import numpy as np
import ast
import sys, os
from scipy.fftpack import fft, ifft
import scipy.io as sio

def loadData(file_path):

    """
    Load a log file. (Use the built-in pickle module)

    Parameters
    ----------
    file_path: str
        path of the log file

    Returns
    -------
    - data: arbitrary object
        arbitrary Python object which contains data
    """
    import pickle
    with open(file_path, 'rb') as _input:
        data = pickle.load(_input)
    return data

def white(w, S0=1.):
    return S0 * np.ones_like(w)
def lorentz(w, S0=1., wc=1., w0=0.):
    """ Lorentzian spectrum. """
    return S0/(2*np.pi*wc)*(1./(1+((w-w0)/wc)**2) + 1./(1+((w+w0)/wc)**2))

def genProcess(times, spectrum, T0):
    """ Generate one realization of a Gaussian random process. """

    # times:    numpy array containing all the times for which to generate the process
    # spectrum: numpy array containing the spectrum at all the harmonics omega_k = 2*pi*k/T0
    # T0:       period of the noise process. Has to be much longer than the duration of an
    #           experimental run
    # returns:  numpy array containing all the values of the random process for the times
    #           specified in input

    # Check that T0 is longer than the largest specified time
    if (T0 < max(times)):
        print("WARNING: T0 should be the longest timescale.")

    # Number of frequencies and times
    nfreqs = max(np.shape(spectrum))
    ntimes = max(np.shape(times))

    # Generation of harmonics and random Fourier coefficients
    vec_omega = 2*np.pi/T0 * np.arange(1,nfreqs+1)
    vec_sigma = np.sqrt(2.*spectrum/T0)
    vec_a = vec_sigma * np.random.normal(size=nfreqs)
    vec_b = vec_sigma * np.random.normal(size=nfreqs)

    # Matrix with times in rows
    mat_times = np.array(np.repeat(times[np.newaxis],nfreqs,0))

    # Matrix with frequencies in columns
    mat_freqs = np.array(np.repeat(np.transpose(vec_omega[np.newaxis]),ntimes,1))

    # Sum up the Fourier series
    mat_cos = np.cos(mat_freqs * mat_times);
    cos_term = np.dot(np.reshape(vec_a,(1,nfreqs)), mat_cos)
    mat_sin = np.sin(mat_freqs * mat_times);
    sin_term = np.dot(np.reshape(vec_b,(1,nfreqs)), mat_sin)

    return (cos_term + sin_term)[0]

class Driver(InstrumentDriver.InstrumentWorker):
    """ This class implements a Single-qubit pulse generator"""

    def performOpen(self, options={}):
        self.vEnvelope = np.array([], dtype=float)
        self.vNoise_Time = np.array([], dtype=float)
        self.vNoise_Time_Modulated = np.array([], dtype=float)
        self.vNoise_Freq = np.array([], dtype=float)
        self.vNoise_Freq_FFT = np.array([], dtype=float)
        self.vTime = np.array([], dtype=float)
        self.vFreq = np.array([], dtype=float)
        self.vFreq_AWG_frame = np.array([], dtype=float)
        self.vBinEdges = np.array([], dtype=float)
        self.vHistogram = np.array([], dtype=float)
        self.PrevNoiseIndex = -9.99e-9
        self.saved_data = None
        self.original_S0 = 0

    def performSetValue(self, quant, value, sweepRate=0.0, options={}):
        """Perform the Set Value instrument operation. This function should
        return the actual value set by the instrument"""
        # do nothing, just return value
        if quant.name == 'Read file':
            self.Readfile()
        elif quant.name == 'Read file when start':
            if (value == True):
                self.Readfile()
        return value

    def Readfile(self):
        # self.log('Read File')
        file_path = self.getValue("File path")
        self.saved_data = loadData(file_path)
        self.setValue("T0", self.saved_data["T0"])
        self.setValue("High Cutoff Freq.", self.saved_data["High Cutoff Freq."])
        self.original_S0 = self.saved_data["S0"]
        self.setValue("Noise Power", self.saved_data["S0"])
        self.setValue("Center Freq.", self.saved_data["Center Freq."])
        self.setValue("HWHM", self.saved_data["HWHM"])
        self.reportStatus('Successfully Read Data File!')

    def performGetValue(self, quant, options={}):
        """Perform the Get Value instrument operation"""
        # check type of quantity
        if quant.isVector():
            noise_indx = self.getValue("Noise Index")
            self.calculateNoise()
            # if self.PrevNoiseIndex != noise_indx:
            #     self.calculateNoise()
            #     self.PrevNoiseIndex = noise_indx

            # traces, check if waveform needs to be re-calculated
            if self.isConfigUpdated():
                self.calculateWaveform()

            # if ((quant.name == 'Trace - Modulated Noise (Time-Domain)') & (len(self.vEnvelope) != len(self.vNoise_Time))):
            #     self.calculateNoise()

            # get correct data and return as trace dict
            vData = self.getWaveformFromMemory(quant)
            length = len(vData)
            if (self.getValue("Turn Off Noise") == True):
                vData = np.zeros(length)
            if quant.name == 'Trace - Noise (Histogram)':
                value = quant.getTraceDict(vData, x=self.vBinEdges)
            elif quant.name == 'Trace - Noise (Freq-Domain, FFT)':
                value = quant.getTraceDict(vData, x=self.vFreq_AWG_frame)
            elif quant.name == 'Trace - Noise (Freq-Domain, Original)':
                value = quant.getTraceDict(vData, x=self.vFreq)
            # elif quant.name == 'Trace - Envelope':
            #     dt = 1/self.getValue('Sample rate')
            #     value = quant.getTraceDict(self.vEnvelope, dt = dt)
            else:
                dt = 1/self.getValue('Sample rate')
                value = quant.getTraceDict(vData, dt=dt)
        else:
            # for all other cases, do nothing
            value = quant.getValue()

        return value

    def getWaveformFromMemory(self, quant):

        dTrace = { 'Trace - Noise (Time-Domain)': self.vNoise_Time,
                  'Trace - Noise (Freq-Domain, Original)': self.vNoise_Freq,
                  # 'Trace - Noise (Time-Domain, RAW)': self.vNoise_Time_RAW,
                  'Trace - Noise (Freq-Domain, FFT)': self.vNoise_Freq_FFT,
                  'Trace - Noise (Histogram)': self.vHistogram,
                  'Trace - Envelope': self.vEnvelope,
                  'Trace - Modulated Noise (Time-Domain)': self.vNoise_Time * self.vEnvelope if (len(self.vNoise_Time) == len(self.vEnvelope)) else []
                  # 'Trace - Noise (Time-Domain, Histogram, RAW)': self.vTime_Histogram_RAW,
                  }
        vData = dTrace[quant.name]
        return vData
    def generateNoise(self):
        """Get the waveform of a given noise"""
        # get params
        nPoints = float(self.getValue('Number of points'))
        dSampleRate = self.getValue('Sample rate')
        str_NoiseType = self.getValue('Noise type')
        index = int(self.getValue('Noise Index'))
        use_interpolation = self.getValue('Use Interpolation')

        if (str_NoiseType == 'Custom'):
            vec_time = np.linspace(0, (nPoints-1)*1/dSampleRate,nPoints)
            vec_freq_AWG_frame = np.linspace(0.0, 1.0/(2.0/dSampleRate), nPoints*0.5)
            if (self.saved_data is None):
                self.reportStatus('No Data Available')
                return
            vec_T_noise = np.zeros_like(vec_time)
            index = np.clip(index, 0, len(self.saved_data['mat_T_noise'])-1)
            self.log('shape of mat_T_noise: ' + str(self.saved_data['mat_T_noise'].shape))
            vec_T_noise_data = np.copy(self.saved_data['mat_T_noise'][index]) * np.sqrt(self.getValue('Noise Power')/self.original_S0)
            self.log( self.saved_data['mat_T_noise'][index])
            if (len(vec_T_noise_data) >= len(vec_time)):
                vec_T_noise = vec_T_noise_data[:len(vec_time)]
            else:
                N_copy = int(len(vec_time) / len(vec_T_noise_data))
                for i in range(N_copy):
                    vec_T_noise[i*len(vec_T_noise_data):(i+1)*len(vec_T_noise_data)] = vec_T_noise_data
                vec_T_noise[(i+1)*len(vec_T_noise_data):] = vec_T_noise_data[0:len(vec_time)-(i+1)*len(vec_T_noise_data)]
            vec_F_noise = self.saved_data['vec_F_noise']
            vec_freq = self.saved_data['vec_freq']

        else:

            np.random.seed(int(index))
            T0 = self.getValue('T0')
            dHighCutoffFreq = self.getValue('High Cutoff Freq.')
            if use_interpolation:
                interval = self.getValue('Interpolation Interval')
                dSampleRate = int(dSampleRate / interval)
                nPoints = int(nPoints / self.getValue('Sample rate') * dSampleRate)
            # self.log(T0, dSampleRate, nPoints)
            vec_time = np.linspace(0, (nPoints-1)*1/dSampleRate,nPoints)
            vec_freq_AWG_frame = np.linspace(0.0, 1.0/(2.0/dSampleRate), nPoints*0.5)

            S0 = self.getValue('S0')
            freq_step = 1/T0
            N_freq = int(dHighCutoffFreq / freq_step)
            vec_freq = np.linspace(freq_step, dHighCutoffFreq, N_freq)

            if str_NoiseType == 'White':
                """Generate Gaussian Noise Signal(mean = 0, sigma = sqrt(Power Spectral Density))"""
                vec_F_noise = white(2 * np.pi * vec_freq, S0= S0)
                vec_T_noise = genProcess(vec_time, vec_F_noise, T0)

            elif str_NoiseType == 'Squared-Gaussian':
                """Generate Squared Gaussian Noise Signal"""
                vec_F_noise = white(2 * np.pi * vec_freq, S0= S0)
                vec_T_noise = genProcess(vec_time, vec_F_noise, T0) **2

            elif str_NoiseType == 'Lorentzian':
                f0 = self.getValue('Center Freq.')
                fwidth = self.getValue('HWHM')
                vec_F_noise = lorentz(2 * np.pi * vec_freq , S0= S0, wc = 2 * np.pi * fwidth , w0 = 2 * np.pi * f0)
                vec_T_noise = genProcess(vec_time, vec_F_noise, T0)

        #Filter noise in Time-Window
        dNoiseStart = self.getValue('Noise Start Time')
        dNoiseEnd = self.getValue('Noise End Time')
        index_Empty = np.where((vec_time < dNoiseStart) | (vec_time > dNoiseEnd))
        vec_T_noise[index_Empty] = 0

        bln_hist_FFT = self.getValue('Generate Histogram and FFT')
        d_bincounts = self.getValue('Histogram Bin Counts')
        vec_F_noise_FFT = []
        vec_histogram = []
        vec_binedges = []
        if (use_interpolation == True):
            nPoints = float(self.getValue('Number of points'))
            dSampleRate = self.getValue('Sample rate')
            vec_time_interp = np.linspace(0, (nPoints-1)*1/dSampleRate,nPoints)
            vec_freq_AWG_frame = np.linspace(0.0, 1.0/(2.0/dSampleRate), nPoints*0.5)
            interp_T_noise = np.interp(vec_time_interp, vec_time, vec_T_noise)
            vec_time = np.copy(vec_time_interp)
            vec_T_noise = np.copy(interp_T_noise)

        if (bln_hist_FFT):
            vec_histogram, vec_binedges = np.histogram(vec_T_noise, bins = int(d_bincounts))
            vec_F_noise_FFT = 2.0 / len(vec_time) * np.abs(fft(vec_T_noise).flatten())

        else:
            vec_F_noise_FFT = np.zeros_like(vec_freq_AWG_frame)
        return ({
                 "vec_time": vec_time,
                 "vec_freq": vec_freq,
                 "vec_freq_AWG_frame": vec_freq_AWG_frame,
                 "vec_T_noise": vec_T_noise,
                 "vec_F_noise_FFT": vec_F_noise_FFT[:len(vec_time)//2],
                 "vec_F_noise":vec_F_noise,
                 "vec_histogram": vec_histogram,
                 "vec_binedges": vec_binedges[0:int(d_bincounts)],
                 })

    def calculateNoise(self):
        """Generate noise waveform"""
        # get config values
        nPoints = int(self.getValue('Number of points'))
        sampleRate = self.getValue('Sample rate')
        # start with allocating time and amplitude vectors
        self.vTime = np.arange(nPoints, dtype=float)/sampleRate

        # get noise waveform
        result= self.generateNoise()
        if (result is None):
            return
        self.vTime = result["vec_time"]
        self.vFreq = result["vec_freq"]
        self.vFreq_AWG_frame = result["vec_freq_AWG_frame"]
        self.vBinEdges = result["vec_binedges"]
        self.vNoise_Time = result["vec_T_noise"]
        self.vNoise_Freq = result["vec_F_noise"]
        self.vNoise_Freq_FFT = result["vec_F_noise_FFT"]
        self.vHistogram = result["vec_histogram"]

    def calculateWaveform(self):
        """Generate waveforms, including pre-pulses, readout and gates"""
        # get config values
        nPoints = int(self.getValue('Number of points'))
        sampleRate = self.getValue('Sample rate')
        firstDelay = self.getValue('First pulse delay')
        # start with allocating time and amplitude vectors
        self.vTime = np.arange(nPoints, dtype=float)/sampleRate
        # create list of output vectors
        self.vEnvelope = np.zeros_like(self.vTime)
        # go on depending on with sequence
        self.generateSequence(startTime=firstDelay)

    def generateSequence(self, startTime):
        # get config values
        sSequence = self.getValue('Sequence')
        nPulses = int(self.getValue('# of pulses'))
        seqPeriod = self.getValue('Pulse period')
        # go on depending on waveform
        if sSequence == 'CP/CPMG':
            # get length of actual pulses
            dPulseT1 = self.getPulseDuration(1)
            dPulseT2 = self.getPulseDuration(2)
            dPulseTot = 2*dPulseT1 + dPulseT2*nPulses
            # add the first pi/2 pulses
            self.addPulse(1, startTime + dPulseT1/2)
            # add more pulses
            if nPulses <= 0:
                # no pulses = ramsey
                vTimePi = []
                # second pi/2 pulse
                self.addPulse(1, startTime + seqPeriod + dPulseTot - dPulseT1/2)
            elif nPulses == 1:
                # one pulse, echo experiment
                vTimePi = [startTime + dPulseT1 + seqPeriod/2 + dPulseT2/2]
                # second pi/2 pulse
                self.addPulse(1, startTime + seqPeriod + dPulseTot - dPulseT1/2)
            elif nPulses > 1:
                # figure out timing of pi pulses
                vTimePi = startTime  + dPulseT1 + seqPeriod/2 + dPulseT2/2 + \
                          (seqPeriod + dPulseT2)*np.arange(nPulses)
                # second pi/2 pulse
                self.addPulse(1, startTime + nPulses*seqPeriod + dPulseTot - dPulseT1/2)
            # add pi pulses, one by one
            for dTimePi in vTimePi:
                self.addPulse(2, dTimePi)

        elif sSequence == 'Generic sequence':
            # generic pulse sequence, add the pulses specified in the pulse list
            t = startTime
            for n in range(nPulses):
                pulseType = 1 + (n % 8)
                # get length of current pulse
                dPulseT = self.getPulseDuration(pulseType)
                self.addPulse(pulseType, t + dPulseT/2)
                # add spacing as defined for this pulse
                t += dPulseT + self.getValue('Spacing #%d' % (pulseType))

    def addPulse(self, nType, dTime, nOutput=None, bTimeStart=False, phase=None):
        """Add pulse to waveform"""
        vTime, vPulse, vIndx = self.getPulseEnvelope(nType, dTime, bTimeStart)
        if len(vTime) == 0:
            return
        self.vEnvelope[vIndx] += vPulse

    def getPulseEnvelope(self, nType, dTime, bTimeStart=False):
        """Get pulse envelope for a given pulse"""
        sPulseType = self.getValue('Pulse type')
        dSampleRate = self.getValue('Sample rate')
        truncRange = self.getValue('Truncation range')
        start_at_zero = self.getValue('Start at zero')
        # get pulse params
        dAmp = self.getValue('Amplitude #%d' % nType)
        dWidth = self.getValue('Width #%d' % nType)
        dPlateau = self.getValue('Plateau #%d' % nType)
        # get pulse width
        if sPulseType == 'Square':
            dTotTime = dWidth+dPlateau
        elif sPulseType == 'Ramp':
            dTotTime = 2*dWidth+dPlateau
        elif sPulseType == 'Gaussian':
            dTotTime = truncRange*dWidth + dPlateau
        # shift time to mid point if user gave start point
        if bTimeStart:
            dTime = dTime + dTotTime/2
        # get the range of indices in use
        vIndx = np.arange(max(np.round((dTime-dTotTime/2)*dSampleRate), 0),
                          min(np.round((dTime+dTotTime/2)*dSampleRate), len(self.vTime)))
        vIndx = np.int0(vIndx)
        self.log(len(vIndx), len(self.vTime))
        # calculate time values for the pulse indices
        vTime = vIndx/dSampleRate
        # calculate the actual value for the selected indices
        if sPulseType == 'Square':
            vPulse = (vTime >= (dTime-(dWidth+dPlateau)/2)) & \
                 (vTime < (dTime+(dWidth+dPlateau)/2))
        elif sPulseType == 'Ramp':
            # rising and falling slopes
            vRise = (vTime-(dTime-dPlateau/2-dWidth))/dWidth
            vRise[vRise<0.0] = 0.0
            vRise[vRise>1.0] = 1.0
            vFall = ((dTime+dPlateau/2+dWidth)-vTime)/dWidth
            vFall[vFall<0.0] = 0.0
            vFall[vFall>1.0] = 1.0
            vPulse = vRise * vFall
#            vPulse = np.min(1, np.max(0, (vTime-(dTime-dPlateau/2-dWidth))/dWidth)) * \
#               np.min(1, np.max(0, ((dTime+dPlateau/2+dWidth)-vTime)/dWidth))
        elif sPulseType == 'Gaussian':
            # width is two times std
            #dStd = dWidth/2;
            # alternate def; std is set to give total pulse area same as a square
            dStd = dWidth/np.sqrt(2*np.pi)
            # cut the tail part and increase the amplitude, if necessary
            dOffset = 0
            if dPlateau > 0:
                # add plateau
                vPulse = (vTime >= (dTime-dPlateau/2)) & \
                    (vTime < (dTime+dPlateau/2))
                if dStd > 0:
                    # before plateau
                    vPulse = vPulse + (vTime < (dTime-dPlateau/2)) * \
                        (np.exp(-(vTime-(dTime-dPlateau/2))**2/(2*dStd**2))-dOffset)/(1-dOffset)
                    # after plateau
                    vPulse = vPulse + (vTime >= (dTime+dPlateau/2)) * \
                        (np.exp(-(vTime-(dTime+dPlateau/2))**2/(2*dStd**2))-dOffset)/(1-dOffset)
            else:
                if dStd > 0:
                    vPulse = (np.exp(-(vTime-dTime)**2/(2*dStd**2))-dOffset)/(1-dOffset)
                else:
                    vPulse = np.zeros_like(vTime)
#        # add the pulse to the previous ones
#        vY[iPulse] = vY[iPulse] + dAmp * vPulse
        vPulse = dAmp * vPulse
        if start_at_zero:
            vPulse = vPulse - vPulse.min()
            vPulse = vPulse/vPulse.max()*dAmp
        # return both time, envelope, and indices
        return (vTime, vPulse, vIndx)

    def getPulseDuration(self, nType):
        """Get total pulse duration waveform, for timimg purposes"""
        # check if edge-to-edge
        if self.getValue('Edge-to-edge pulses'):
            width = self.getValue('Width #%d' % nType)
            plateau = self.getValue('Plateau #%d' % nType)
            pulseEdgeWidth = self.getValue('Edge position')
            return pulseEdgeWidth * width + plateau
        else:
            return 0.0

max_stars_repo_name: roniwinik/Drivers

max_issues_repo_name: roniwinik/Drivers

max_forks_repo_name: roniwinik/Drivers


text:
import os
import numpy as np
import PIL
import random

from torch.utils.data import Dataset
from torchvision.datasets import CIFAR100
from torchvision.transforms import Compose, Lambda, ToTensor, RandomHorizontalFlip, Pad, RandomCrop, Normalize
from PIL import Image

class CustomCIFAR100(Dataset):
    def __init__(self, train=True):
        super(CustomCIFAR100, self).__init__()
        self.cifar_100 = CIFAR100(root='datasets', train = train, download=True)
        #  self._cifar10_train => index[0] = data, index[1] = label, format = PIL

        tensors = list()
        for i in range(len(self.cifar_100)):
            tensors.append(ToTensor()(self.cifar_100[i][0]).numpy())  # cifar100
            # ToTensor 는 들어온 이미지나 ndarray(H, W, C)를 (C, H, W)의 텐서로 바꿔서 돌려줌
        mean = np.mean(tensors, axis=(0,2,3))  # 채널 별 평균 구하기. ToTensor를 썼기 때문에 index=1이 채널
        std = np.std(tensors, axis=(0,2,3))  # 채널 별 표준편차 구하기. ToTensor를 썼기 때문에 index=1이 채널
        print("mean: {}, std: {}".format(mean, std))

        transform = [RandomHorizontalFlip()]
        transform += [Pad(4), RandomCrop(32)]
        transform += [ToTensor(), Normalize(mean=mean, std=std)]
        self.transform = Compose(transform)

    def __getitem__(self, index):  # __getitem__ : index로 value를 얻을 수 있게 하는 super method
        tensor, label = self.transform(self.cifar_100[index][0]), self.cifar_100[index][1]
        return tensor, label

    def __len__(self):  # __getitem__ : len()로 데이터의 길이를 얻을 수 있게 하는 super method
        return len(self.cifar_100)

if __name__=='__main__':
    import os
    dataset = CustomCIFAR100()
    print(dataset[0])

max_stars_repo_name: KWYi/PytorchStudy

max_issues_repo_name: KWYi/PytorchStudy

max_forks_repo_name: KWYi/PytorchStudy


text:
#!/usr/bin/env python
"""Infers correlations between direct effects and population effects, and between direct effects and average non-transmitted coefficients (NTCs).

Minimally: the script requires summary statistics as output by snipar's gwas.py script, and
either LD-scores (as output by snipar's ibd.py script or LDSC) or .bed files from
which LD-scores can be computed

Args:
@parser@

Results:
    correlations
        A text file containing the estimated correlations and their standard errors.

"""
import numpy as np
import argparse
from snipar.correlate import *
from numba import set_num_threads
from numba import config as numba_config
from snipar.utilities import *
from snipar.utilities import get_parser_doc
parser = argparse.ArgumentParser()
parser.add_argument('sumstats', type=str, help='Address of sumstats files in SNIPar sumstats.gz text format (without .sumstats.gz suffix). If there is a @ in the address, @ is replaced by the chromosome numbers in chr_range (optional argument)')
parser.add_argument('--chr_range',
                        type=parseNumRange,
                        nargs='*',
                        action=NumRangeAction,
                        help='number of the chromosomes to be imputed. Should be a series of ranges with x-y format or integers.', default=None)
parser.add_argument('out',type=str,help='Prefix for output file(s)')
parser.add_argument('--ldscores',type=str,help='Address of ldscores as output by LDSC',default=None)
parser.add_argument('--bed', type=str,
                    help='Address of observed genotype files in .bed format (without .bed suffix). If there is a # in the address, # is replaced by the chromosome numbers in the range of 1-22.',
                    default=None)
parser.add_argument('--threads',type=int,help='Number of threads to use for IBD inference. Uses all available by default.',default=None)
parser.add_argument('--min_maf',type=float,help='Ignore SNPs with minor allele frequency below min_maf (default 0.05)', default=0.05)
parser.add_argument('--corr_filter',type=float,help='Filter out SNPs with outlying sampling correlations more than corr_filter SDs from mean (default 6)',default=6.0)
parser.add_argument('--n_blocks',type=int,help='Number of blocks to use for block-jacknife variance estimate (default 200)',default=200)
parser.add_argument('--save_delete',action='store_true',help='Save jacknife delete values',default=False)
parser.add_argument('--ld_wind',type=float,help='The window, in cM, within which LD scores are computed (default 1cM)',default=1.0)
parser.add_argument('--ld_out',type=str,help='Output LD scores in LDSC format to this address',default=None)
__doc__ = __doc__.replace("@parser@", get_parser_doc(parser))
def main(args):
    """"Calling this function with args is equivalent to running this script from commandline with the same arguments.
    Args:
        args: list
            list of all the desired options and arguments. The possible values are all the values you can pass this script from commandline.
    """
    # Set number of threads
    if args.threads is not None:
        num_threads = min([args.threads, numba_config.NUMBA_NUM_THREADS])
    else:
        num_threads = numba_config.NUMBA_NUM_THREADS
    set_num_threads(num_threads)
    print('Number of threads: '+str(num_threads))

    if args.ldscores is None and args.bed is None:
        raise(ValueError('Must provide either LD scores or genotypes to compute ld scores'))
    if args.ldscores is not None and args.bed is not None:
        raise(ValueError('Both LD scores and genotypes provided. Provided LD scores will be used'))

    # Find sumstats files
    sumstats_files, chroms = parse_obsfiles(args.sumstats, obsformat='sumstats.gz', chromosomes=args.chr_range)
    # Read sumstats
    s = read_sumstats_files(sumstats_files, chroms)

    # Filter
    print('Filtering on missing values')
    s.filter_NAs()
    print('Filtering out SNPs with MAF<'+str(args.min_maf))
    s.filter_maf(args.min_maf)
    print('Filtering out SNPs with sampling correlations more than '+str(args.corr_filter)+' SDs from mean')
    s.filter_corrs(args.corr_filter)

    # Get LD scores
    if args.ldscores is not None:
        ld_files, chroms = parse_obsfiles(args.ldscores, obsformat='l2.ldscore.gz', chromosomes=[x for x in range(1,23)])
        s.scores_from_ldsc(ld_files)
        s.filter_NAs()
    elif args.bed is not None:
        bedfiles, chroms = parse_obsfiles(args.bed, obsformat='bed', chromosomes=chroms)
        s.compute_ld_scores(bedfiles, chroms, args.ld_wind, args.ld_out)
        s.filter_NAs()

    # Compute correlations
    print('Using '+str(s.sid.shape[0])+' SNPs to compute correlations')
    r_dir_pop, r_dir_pop_SE, r_dir_pop_delete = s.cor_direct_pop(args.n_blocks)
    print('Correlation between direct and population effects: '+str(round(r_dir_pop,4))+' (S.E. '+str(round(r_dir_pop_SE,4))+')')
    r_dir_avg_NTC, r_dir_avg_NTC_SE, r_dir_avg_NTC_delete = s.cor_direct_avg_NTC(args.n_blocks)
    print('Correlation between direct and average NTCs: '+str(round(r_dir_avg_NTC,4))+' (S.E. '+str(round(r_dir_avg_NTC_SE,4))+')')

    outfile = str(args.out)+'_corrs.txt'
    print('Saving correlation estimates to '+outfile)
    first_col = np.array(['r_direct_population','r_direct_avg_NTC']).reshape((2,1))
    header = np.array(['correlation','est','SE']).reshape((1,3))
    cor_out = np.zeros((2,2))
    cor_out[0,:] = np.array([r_dir_pop,r_dir_pop_SE])
    cor_out[1,:] = np.array([r_dir_avg_NTC,r_dir_avg_NTC_SE])
    outarray = np.vstack((header,np.hstack((first_col,cor_out))))
    np.savetxt(outfile,outarray,fmt='%s')

    if args.save_delete:
        delete_outfile = args.out+'_delete.txt'
        print('Saving jacknife delete values to '+str(delete_outfile))
        delete_out = np.vstack((r_dir_pop_delete,r_dir_avg_NTC_delete)).T
        delete_out = np.vstack((np.array(['r_direct_population','r_direct_avg_NTC']).reshape((1,2)),delete_out))
        np.savetxt(delete_outfile, delete_out,fmt='%s')

if __name__ == "__main__":
    args=parser.parse_args()
    main(args)

max_stars_repo_name: AlexTISYoung/snipar

max_issues_repo_name: AlexTISYoung/snipar

max_forks_repo_name: AlexTISYoung/snipar

reason: import numpy,from numba


text:
import argparse
import os
import shutil
import time
from pathlib import Path

import numpy as np
import cv2
import torch
import torch.backends.cudnn as cudnn

from models.experimental import attempt_load
from utils.datasets import LoadStreams, LoadImages
from utils.general import (
    check_img_size, non_max_suppression, apply_classifier, scale_coords,
    xyxy2xywh, plot_one_box, strip_optimizer, set_logging)
from utils.torch_utils import select_device, load_classifier, time_synchronized

model = 0
device = 0

def mask2gt(mask):
    # instances are encoded as different colors
    obj_ids = np.unique(mask)

    # first id is the background, so remove it
    obj_ids = obj_ids[1:]

    # split the color-encoded mask into a set
    # of binary masks
    masks = mask == obj_ids[:, None, None]

    boxes = []
    # if no objects exists
    if len(obj_ids) == 0:
        return boxes

    for i, obj in enumerate(obj_ids):
        pos = np.where(masks[i])
        xmin, xmax = np.min(pos[1]), np.max(pos[1])
        ymin, ymax = np.min(pos[0]), np.max(pos[0])
        boxes.append([float(xmin)/mask.shape[1], float(ymin)/mask.shape[0],
                      float(xmax)/mask.shape[1], float(ymax)/mask.shape[0]])
    return boxes

def label_unique(mask):
    new_mask = np.zeros(mask.shape)
    contours, _ = cv2.findContours(mask, cv2.RETR_EXTERNAL,cv2.CHAIN_APPROX_NONE)
    for i, cont in enumerate(contours):
        cv2.drawContours(new_mask, [cont], contourIdx=-1, color=(i+1), thickness=-1)
    return new_mask

def detect(batch, offsets):

    img = torch.from_numpy(batch).to(device)
    img = img.half() if half else img.float()  # uint8 to fp16/32
    #img /= 255.0  # 0 - 255 to 0.0 - 1.0
    if img.ndimension() == 3:
        img = img.unsqueeze(0)

    # Inference
    t1 = time_synchronized()
    pred = model(img, augment=opt.augment)[0]

    # Apply NMS
    pred = non_max_suppression(pred, opt.conf_thres, opt.iou_thres, classes=opt.classes, agnostic=opt.agnostic_nms)
    t2 = time_synchronized()

    confs, bbs = [], []
    # Process detections
    for i, det in enumerate(pred):  # detections per image
        if det is not None and len(det):
            off = offsets[i]
            det[:, :4] = scale_coords(img.shape[2:], det[:, :4], (1024,1024,3)).round()

            for d in det:
                x, y, x_, y_, c = int(d[0]), int(d[1]), int(d[2]), int(d[3]), float(d[4])
                bbs.append([x+off[0], y+off[1], x_-x, y_-y])
                confs.append(x)
    return confs, bbs

def data_loader(crop_size = 1024, root_dir = '/home/markpp/datasets/bo/test/', version='mixed', make_gt=False):

    n_width, n_height = (3300+crop_size)//crop_size, (6050+crop_size)//crop_size
    batch_size = n_width * n_height
    print(batch_size)

    items = sorted([f for f in os.listdir(root_dir) if 'item' in f])
    # iterate over each item
    for item in items[:]:
        img_paths = sorted([f for f in os.listdir(os.path.join(root_dir,item,"rgb")) if f.endswith('.jpg')])
        print("{} has {} images".format(item,len(img_paths)))
        # iterate over each image belonging to the given item
        if not os.path.exists(os.path.join(root_dir,item,"results_{}".format(version))):
            os.makedirs(os.path.join(root_dir,item,"results_{}".format(version)))
        if not os.path.exists(os.path.join(root_dir,item,"pred_{}".format(version))):
            os.makedirs(os.path.join(root_dir,item,"pred_{}".format(version)))

        for i,img_file in enumerate(img_paths[:]):
            start = time.time()

            pred_bb, pred_cls, pred_conf = [], [], []
            gt_bb, gt_cls = [], []

            img = cv2.imread(os.path.join(root_dir,item,"rgb",img_file))
            img = cv2.copyMakeBorder(img, top=0, bottom=crop_size, left=0, right=crop_size, borderType=cv2.BORDER_CONSTANT, value=(0,0,0))

            #img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)
            if img is None:
                print("image {} not loaded".format(img_file))
                continue

            if make_gt:
                mask_file = img_file.replace('jpg','png')
                mask = cv2.imread(os.path.join(root_dir,item, "masks",mask_file), cv2.IMREAD_GRAYSCALE)
                if mask is None:
                    mask = np.zeros((img.shape[0],img.shape[1]), dtype=np.uint8)
                    print("empty mask created in place of missing file {}".format(mask_file))
                else:
                    mask = cv2.copyMakeBorder(mask, top=0, bottom=crop_size, left=0, right=crop_size, borderType=cv2.BORDER_CONSTANT, value=0)

                new_mask  = label_unique(mask)
                #cv2.imwrite(os.path.join(root_dir,item,"results",mask_file),new_mask)

                gts = mask2gt(new_mask)
                for gt in gts:
                    gt_bb.append(gt) #xmin, ymin, xmax, ymax
                    gt_cls.append(0)

            offsets = []
            batch = np.empty((batch_size,3,crop_size,crop_size))
            idx = 0
            # divide image into appropritate crops and detect defects
            for x in range(0,img.shape[1]-crop_size,crop_size):
                for y in range(0,img.shape[0]-crop_size,crop_size):
                    crop = img[y:y+crop_size, x:x+crop_size]
                    batch[idx] = np.rollaxis(crop, 2, 0) / 255
                    idx += 1
                    offsets.append([x,y])
                    if idx == batch_size:
                        #detect(batch, offsets)

                        confs, bbs = detect(batch, offsets)

                        for conf, bb in zip(confs, bbs):
                            left, top, width, height = bb[0], bb[1], bb[2], bb[3] # add crop offset
                            img = cv2.rectangle(img, (left, top), (left + width, top + height), (0,255,0), 2)
                            cv2.putText(img, "{:.2f}".format(conf), (left - 10, top - 10), cv2.FONT_HERSHEY_SIMPLEX, 0.5, (0,0,255), 2)

                            pred_bb.append([float(left)/img.shape[1], float(top)/img.shape[0],
                                            float(left+width)/img.shape[1], float(top+height)/img.shape[0]])
                            pred_cls.append(0)
                            pred_conf.append(conf)
                        batch = []
                        offsets = []

            cv2.imwrite(os.path.join(root_dir,item,"results_{}".format(version),img_file),img)

            pred_name = img_file.replace('jpg','txt')
            pred_path = os.path.join(root_dir,item,"pred_{}".format(version),pred_name) #os.path.join(root_dir,item,"pred",pred_name)
            if len(pred_bb):
                with open(pred_path, 'w') as pred_file:
                    for conf, bb in zip(pred_conf,pred_bb):
                        pred_file.write('defect {:.3f} {} {} {} {}\n'.format(conf,
                                                                         int(bb[0]*img.shape[1]),
                                                                         int(bb[1]*img.shape[0]),
                                                                         int(bb[2]*img.shape[1]),
                                                                         int(bb[3]*img.shape[0])))
            else:
                np.savetxt(pred_path, [], delimiter=",", fmt='%u')

            if make_gt:
                gt_name = img_file.replace('jpg','txt')
                gt_path = os.path.join(root_dir,item,"gt",gt_name) #os.path.join(root_dir,item,"gt",gt_name)
                if len(gt_bb):
                    with open(gt_path, 'w') as gt_file:
                        for bb in gt_bb:
                            gt_file.write('defect {} {} {} {}\n'.format(int(bb[0]*img.shape[1]),
                                                                        int(bb[1]*img.shape[0]),
                                                                        int(bb[2]*img.shape[1]),
                                                                        int(bb[3]*img.shape[0])))
                else:
                    np.savetxt(gt_path, [], delimiter=",", fmt='%u')

            end = time.time()
            print("{} took {}".format(img_file, end - start))

#python3 my_detect.py --weights best.pt --img-size 1024

if __name__ == '__main__':
    parser = argparse.ArgumentParser()
    parser.add_argument('--weights', nargs='+', type=str, default='yolov5s.pt', help='model.pt path(s)')
    parser.add_argument('--source', type=str, default='inference/images', help='source')  # file/folder, 0 for webcam
    parser.add_argument('--img-size', type=int, default=640, help='inference size (pixels)')
    parser.add_argument('--conf-thres', type=float, default=0.1, help='object confidence threshold')
    parser.add_argument('--iou-thres', type=float, default=0.25, help='IOU threshold for NMS')
    parser.add_argument('--device', default='', help='cuda device, i.e. 0 or 0,1,2,3 or cpu')
    parser.add_argument('--view-img', action='store_true', help='display results')
    parser.add_argument('--save-txt', action='store_true', help='save results to *.txt')
    parser.add_argument('--save-conf', action='store_true', help='save confidences in --save-txt labels')
    parser.add_argument('--save-dir', type=str, default='inference/output', help='directory to save results')
    parser.add_argument('--classes', nargs='+', type=int, help='filter by class: --class 0, or --class 0 2 3')
    parser.add_argument('--agnostic-nms', action='store_true', help='class-agnostic NMS')
    parser.add_argument('--augment', action='store_true', help='augmented inference')
    parser.add_argument('--update', action='store_true', help='update all models')
    opt = parser.parse_args()
    print(opt)

    with torch.no_grad():

        device = select_device(opt.device)
        half = device.type != 'cpu'  # half precision only supported on CUDA

        # Load model
        model = attempt_load(opt.weights, map_location=device)  # load FP32 model
        imgsz = check_img_size(opt.img_size, s=model.stride.max())  # check img_size
        if half:
            model.half()  # to FP16

        if opt.update:  # update all models (to fix SourceChangeWarning)
            for opt.weights in ['yolov5s.pt', 'yolov5m.pt', 'yolov5l.pt', 'yolov5x.pt']:
                data_loader(crop_size = opt.img_size)
                #detect()
                strip_optimizer(opt.weights)
        else:
            data_loader(crop_size = opt.img_size)
            #detect()


text:
import os
import sys
import h5py
import argparse
import numpy as np
from tqdm import tqdm
opj = os.path.join

from sixd import SixdToolkit

LINEMOD = '/home/penggao/data/sixd/hinterstoisser/test/02'
KPDROOT = '/home/penggao/projects/pose/kp6d/keypoint/data/occ'
NAMES = ('ape', 'bvise', 'bowl', 'camera', 'can', 'cat', 'cup',
                'driller', 'duck', 'eggbox', 'glue', 'holepuncher', 'iron', 'lamp', 'phone')
OCCNAMES = ('ape', 'can', 'cat', 'driller', 'duck',
            'eggbox', 'glue', 'holepuncher')

def parse_arg():
    parser = argparse.ArgumentParser(description='Synthetic data generator')
    parser.add_argument('--seq', type=str, required=True)
    parser.add_argument('--kpnum', choices=[17, 9],
                        type=int, help="Number of keypoints")
    parser.add_argument('--kptype', choices=['sift', 'cluster', 'corner'],
                        type=str, help="Type of keypoints")
    return parser.parse_args()

if __name__ == '__main__':
    args = parse_arg()
    print("[LOG] Preparing h5 for KPD training")
    print("[LOG] Number of keypoints: %d" % args.kpnum)
    print("[LOG] Type of keypoints: %s" % args.kptype)
    print("[LOG] Sequence: %s %s" % (args.seq, NAMES[int(args.seq)-1]))

    if NAMES[int(args.seq) - 1] not in OCCNAMES:
        print("[ERROR] Sequence not exist in occlusion linemod")
        sys.exit()

    bench = SixdToolkit(dataset='hinterstoisser', kpnum=args.kpnum,
                        kptype=args.kptype, is_train=False)

    KPDROOT = opj(KPDROOT, str(args.kpnum), args.kptype, args.seq)

    if os.path.exists(KPDROOT):
        print("[WARNING] Overwriting existing images in %s" % KPDROOT)
        print("[WARNING] Proceed (y/[n])?", end=' ')
        choice = input()
        if choice == 'y':
            os.system('rm -r %s' % KPDROOT)
        else:
            print("[LOG] Terminated by user, exit")
            sys.exit()

    os.makedirs(KPDROOT)
    os.makedirs(opj(KPDROOT, 'train'))
    os.makedirs(opj(KPDROOT, 'eval'))

    print("[LOG] Preparing data")
    with open(opj(LINEMOD, 'train.txt'), 'r') as f:
        trainlists = f.readlines()
    trainlists = [x.strip() for x in trainlists]

    count = 0
    bboxes = {'train': [], 'eval': []}
    kps = {'train': [], 'eval': []}
    imgnames = {'train': [], 'eval': []}
    tbar = tqdm(bench.frames['02'])
    for idx, f in enumerate(tbar):
        tbar.set_description('%04d' % idx)
        tag = 'train' if '%04d' % idx in trainlists else 'eval'

        # annotation
        try:
            annot = f['annots'][f['obj_ids'].index(int(args.seq))]
        except:
            print("Sequence %s not exist in frame %d" % (args.seq, idx))
            continue
        count += 1
        name_chars = []
        for char in '%012d' % int(idx):
            name_chars.append(ord(char))
        bboxes[tag].append(annot['bbox'])
        kps[tag].append(annot['kps'])
        imgnames[tag].append(np.asarray(name_chars))

        # images
        srcimg = f['path']
        dstimg = opj(KPDROOT, tag, '%012d.png' % int(idx))
        os.symlink(srcimg, dstimg)

    for tag in ('train', 'eval'):
        with h5py.File(opj(KPDROOT, 'annot_%s.h5' % tag), "w") as f:
            f.create_dataset("imgname", data=np.vstack(imgnames[tag]))
            f.create_dataset("bndbox", data=np.vstack(
                bboxes[tag]).reshape(-1, 1, 4))
            f.create_dataset("part", data=np.vstack(
                kps[tag]).reshape(-1, args.kpnum, 2))

    print("[LOG] count =", count)
    print("[LOG] Done. H5 file has been generated in %s" % KPDROOT)