text:
# coding=utf-8
# Copyright 2018 The TF-Agents Authors.
#
# 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.

"""Unit tests for PyUniformReplayBuffer and PyHashedReplayBuffer."""

from __future__ import division
from __future__ import unicode_literals

import os

from absl.testing import parameterized
import numpy as np
import tensorflow as tf

from tf_agents.environments import time_step as ts
from tf_agents.environments import trajectory
from tf_agents.policies import policy_step
from tf_agents.replay_buffers import py_hashed_replay_buffer
from tf_agents.replay_buffers import py_uniform_replay_buffer
from tf_agents.specs import array_spec
from tf_agents.utils import nest_utils

def next_dataset_element(test_case, dataset):
  """Utility function to iterate over tf.data.Datasets in both TF 1.x and 2.x.

  TensorFlow 1.x and 2.x have different mechanisms for iterating over elements
  of a tf.data.Dataset. TensorFlow 1.x would require something like:

  itr = tf.data.Dataset.range(10).make_one_shot_iterator()
  get_next = itr.get_next()
  with tf.Session() as sess:
    for _ in range(10):
      item = sess.run(get_next)
      process(item)

  While TensorFlow 2.x enables something simpler like:

  for item in tf.data.Dataset.range(10):
    process(item)

  That simpler latter form is also available in TensorFlow 1.x when running
  with eager execution enabled.

  This function accomodates for the differing styles using:

  next_element = next_dataset_element(self, tf.data.Dataset.range(10))
  for _ in range 10:
    process(next_element())

  Args:
    test_case: The tf.test.TestCase object of the test calling this function.
    dataset: A tf.data.Dataset object.

  Returns:
    A Python function that returns successive elements from dataset on each call
    (using test_case.evaluate() in TensorFlow 1.x).
  """
  if tf.executing_eagerly():
    itr = iter(dataset)
    return lambda: next(itr)
  get_next = tf.compat.v1.data.make_one_shot_iterator(dataset).get_next()
  return lambda: test_case.evaluate(get_next)

class FrameBufferTest(tf.test.TestCase):

  def testFrameBuffer(self):
    fb = py_hashed_replay_buffer.FrameBuffer()
    a = np.random.randint(low=0, high=256, size=[84, 84, 1], dtype=np.uint8)
    h = fb.add_frame(a)
    fb.on_delete([h])
    self.assertEqual(0, len(fb))

    a = np.random.randint(low=0, high=256, size=[84, 84, 1], dtype=np.uint8)
    b = np.random.randint(low=0, high=256, size=[84, 84, 1], dtype=np.uint8)
    fb.add_frame(a)
    h = fb.add_frame(b)
    fb.on_delete([h])
    self.assertEqual(1, len(fb))

class PyUniformReplayBufferTest(parameterized.TestCase, tf.test.TestCase):

  def _create_replay_buffer(self, rb_cls):
    self._stack_count = 4
    self._single_shape = (15, 15, 1)
    shape = (15, 15, self._stack_count)
    observation_spec = array_spec.ArraySpec(shape, np.int32, 'obs')
    time_step_spec = ts.time_step_spec(observation_spec)
    action_spec = policy_step.PolicyStep(array_spec.BoundedArraySpec(
        shape=(), dtype=np.int32, minimum=0, maximum=1, name='action'))
    self._trajectory_spec = trajectory.from_transition(
        time_step_spec, action_spec, time_step_spec)

    self._capacity = 32
    self._replay_buffer = rb_cls(
        data_spec=self._trajectory_spec, capacity=self._capacity)

  def _fill_replay_buffer(self):
    # Generate N frames: the value of pixels is the frame index.
    # The observations will be generated by stacking K frames out of those N,
    # generating some redundancies between the observations.
    single_frames = []
    frame_count = 100
    for k in range(frame_count):
      single_frames.append(np.full(self._single_shape, k, dtype=np.int32))

    # Add stack of frames to the replay buffer.
    time_steps = []
    for k in range(len(single_frames) - self._stack_count + 1):
      observation = np.concatenate(single_frames[k:k + self._stack_count],
                                   axis=-1)
      time_steps.append(ts.transition(observation, reward=0.0))

    self._transition_count = len(time_steps) - 1
    dummy_action = policy_step.PolicyStep(np.int32(0))
    for k in range(self._transition_count):
      self._replay_buffer.add_batch(nest_utils.batch_nested_array(
          trajectory.from_transition(
              time_steps[k], dummy_action, time_steps[k + 1])))

  def _generate_replay_buffer(self, rb_cls):
    self._create_replay_buffer(rb_cls)
    self._fill_replay_buffer()

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testEmptyBuffer(self, rb_cls):
    self._create_replay_buffer(rb_cls=rb_cls)
    ds = self._replay_buffer.as_dataset()
    if tf.executing_eagerly():
      next(iter(ds))
    else:
      get_next = tf.compat.v1.data.make_one_shot_iterator(ds).get_next()
      self.evaluate(get_next)

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testEmptyBufferBatchSize(self, rb_cls):
    self._create_replay_buffer(rb_cls=rb_cls)
    ds = self._replay_buffer.as_dataset(sample_batch_size=2)
    if tf.executing_eagerly():
      next(iter(ds))
    else:
      get_next = tf.compat.v1.data.make_one_shot_iterator(ds).get_next()
      self.evaluate(get_next)

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testEmptyBufferNumSteps(self, rb_cls):
    self._create_replay_buffer(rb_cls=rb_cls)
    ds = self._replay_buffer.as_dataset(num_steps=2)
    if tf.executing_eagerly():
      next(iter(ds))
    else:
      get_next = tf.compat.v1.data.make_one_shot_iterator(ds).get_next()
      self.evaluate(get_next)

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testReplayBufferCircular(self, rb_cls):
    self._generate_replay_buffer(rb_cls=rb_cls)

    # Since data is added in a circular way, we know that frames sampled from
    # the replay buffer should not have values below a given threshold.
    ds = self._replay_buffer.as_dataset()
    next_trajectory = next_dataset_element(self, ds)
    min_value = self._transition_count - self._capacity
    for _ in range(200):
      traj = next_trajectory()
      self.assertLessEqual(min_value, traj.observation[0, 0, 0])
      self.assertAllEqual(traj.observation[:, :, 0] + 1,
                          traj.observation[:, :, 1])
      self.assertAllEqual(traj.observation[:, :, 0] + 2,
                          traj.observation[:, :, 2])
      self.assertAllEqual(traj.observation[:, :, 0] + 3,
                          traj.observation[:, :, 3])

  def testSampleDoesNotCrossHead(self):
    np.random.seed(12345)

    data_spec = array_spec.ArraySpec((), np.int32)
    replay_buffer = py_uniform_replay_buffer.PyUniformReplayBuffer(
        data_spec=data_spec, capacity=10)

    # Seed RB with 5 elements to move head to position 5.
    for _ in range(5):
      replay_buffer.add_batch(np.array([0]))

    # Fill RB with elements 0-9.
    for i in range(10):
      replay_buffer.add_batch(np.array([i]))

    # Sample with num_steps = 2. We should never sample (9, 0) since this is an
    # invalid transition. With 1000 samples, the probability of sampling (9, 0)
    # if it were not protected against would be (1 - (9/10)^10000) ~= 1.
    sample_frequency = [0 for _ in range(10)]
    for _ in range(10000):
      (first, second) = replay_buffer.get_next(num_steps=2, time_stacked=False)
      self.assertNotEqual(np.array(9), first)
      self.assertNotEqual(np.array(0), second)
      sample_frequency[first] += 1

    # 0-9 should all have been sampled about 10000/9 ~= 1111. We allow a delta
    # of 150 off of 1111 -- the chance each sample frequency is within this
    # range is 99.9998% (computed using the pmf of the binomial distribution).
    # And since we fix the random seed, this test is repeatable.
    for i in range(9):
      self.assertAlmostEqual(10000 / 9, sample_frequency[i], delta=150)

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testSampleBatches(self, rb_cls):
    self._generate_replay_buffer(rb_cls=rb_cls)

    ds = self._replay_buffer.as_dataset(sample_batch_size=5)
    next_trajectory = next_dataset_element(self, ds)
    self.assertEqual(list(ds.output_shapes.observation), [5, 15, 15, 4])
    self.assertEqual(list(ds.output_shapes.action), [5])
    traj = next_trajectory()
    self.assertEqual(traj.observation.shape, (5, 15, 15, 4))
    self.assertEqual(traj.step_type.shape, (5,))

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testSampleBatchesWithNumSteps(self, rb_cls):
    self._generate_replay_buffer(rb_cls=rb_cls)

    ds = self._replay_buffer.as_dataset(sample_batch_size=5, num_steps=3)
    self.assertEqual(list(ds.output_shapes.observation), [5, 3, 15, 15, 4])
    self.assertEqual(list(ds.output_shapes.action), [5, 3])
    next_trajectory = next_dataset_element(self, ds)

    traj = next_trajectory()
    self.assertEqual(traj.observation.shape, (5, 3, 15, 15, 4))
    self.assertEqual(traj.action.shape, (5, 3))

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testNumStepsNoBatching(self, rb_cls):
    self._generate_replay_buffer(rb_cls=rb_cls)

    ds = self._replay_buffer.as_dataset(num_steps=3)
    self.assertEqual(list(ds.output_shapes.observation), [3, 15, 15, 4])
    self.assertEqual(list(ds.output_shapes.action), [3])
    next_trajectory = next_dataset_element(self, ds)

    traj = next_trajectory()
    self.assertEqual(traj.observation.shape, (3, 15, 15, 4))
    self.assertEqual(traj.action.shape, (3,))

  @parameterized.named_parameters(
      [('WithoutHashing', py_uniform_replay_buffer.PyUniformReplayBuffer),
       ('WithHashing', py_hashed_replay_buffer.PyHashedReplayBuffer)])
  def testCheckpointable(self, rb_cls):
    self._generate_replay_buffer(rb_cls=rb_cls)
    self.assertEqual(32, self._replay_buffer.size)

    with self.cached_session():
      directory = self.get_temp_dir()
      prefix = os.path.join(directory, 'ckpt')
      saver = tf.train.Checkpoint(rb=self._replay_buffer)
      save_path = saver.save(prefix)

      loaded_rb = (
          rb_cls(data_spec=self._trajectory_spec, capacity=self._capacity))
      loader = tf.train.Checkpoint(rb=loaded_rb)
      loader.restore(save_path).initialize_or_restore()
      self.assertEqual(32, loaded_rb.size)

      # Check that replay buffer contains the same items as before
      ds = loaded_rb.as_dataset()
      next_trajectory = next_dataset_element(self, ds)
      min_value = self._transition_count - self._capacity
      for _ in range(200):
        traj = next_trajectory()
        self.assertLessEqual(min_value, traj.observation[0, 0, 0])
        self.assertAllEqual(traj.observation[:, :, 0] + 1,
                            traj.observation[:, :, 1])
        self.assertAllEqual(traj.observation[:, :, 0] + 2,
                            traj.observation[:, :, 2])
        self.assertAllEqual(traj.observation[:, :, 0] + 3,
                            traj.observation[:, :, 3])

if __name__ == '__main__':
  tf.test.main()


text:
import numpy as np

import astropy.units as u
from astropy.coordinates import ICRS, Galactic, SkyCoord

from pyspi.spi_pointing import SPIPointing
from pyspi.spi_frame import SPIFrame

from pyspi.io.package_data import get_path_of_data_file

def test_spi_pointing_constructor():

    spi_pointing_file = get_path_of_data_file('sc_orbit_param.fits.gz')

    spi_pointing = SPIPointing(spi_pointing_file)

    # assert that the misalignment matrix has not changed

    default_pointing_matrix = np.array([
    [ 9.99997136e-01,  2.65497922e-04, -2.37837002e-03],
    [-2.69012020e-04,  9.99998873e-01, -1.47732581e-03],
    [ 2.37797511e-03,  1.47796139e-03,  9.99996080e-01]
    ])

    assert  np.allclose(spi_pointing._misalignment_matrix,default_pointing_matrix)

# assert the correct shape of the of all the pointing matrices

    assert spi_pointing.sc_matrix.shape == (417,3,3)

    # assert that the test ra,dec point have not changed

    first_sc_points = spi_pointing.sc_points[0]

    # assert  np.allclose(first_sc_points['scx_dec'] ,50.62)
    # assert  np.allclose(first_sc_points['scx_ra'] , 206.32)
    #
    # assert  np.allclose(first_sc_points['scy_dec'] , 32.53)
    # assert np.allclose(first_sc_points['scy_ra'], 65.30, to)
    #
    # assert np.allclose(first_sc_points['scz_dec'],  -19.65)
    # assert  np.allclose(first_sc_points['scz_ra'], 142.13)

def test_spi_frame():

    spi_pointing_file = get_path_of_data_file('sc_orbit_param.fits.gz')

    spi_pointing = SPIPointing(spi_pointing_file)

    # get the first pointing points
    scp = spi_pointing.sc_points[0]

    # construct a SPI Frame

    spi_frame = SPIFrame(lat=0 * u.deg, lon=0 * u.deg, **scp)

    # test that transformation to ICRS make sense
    icrs_along_x_dir = spi_frame.transform_to(ICRS)

    assert np.allclose(icrs_along_x_dir.ra.value, scp['scx_ra'])

    assert np.allclose(icrs_along_x_dir.dec.value, scp['scx_dec'])

    # assert that angular separation is conserved

    source = SkyCoord(ra=197.075, dec=58.98, unit='deg', frame='icrs')
    center = SkyCoord(ra=206.32522477, dec=50.62928435, unit='deg', frame='icrs')

    icrs_sep = source.separation(center).value

    source_spi = source.transform_to(spi_frame)

    center_spi = center.transform_to(spi_frame)

    spi_sep = spi_frame.separation(source_spi).value

    assert np.allclose(icrs_sep, spi_sep)

    # assert that transforms in and out of frame are conserved

ext: py

max_stars_repo_name: mpleinti/pyspi

max_issues_repo_name: mpleinti/pyspi

max_forks_repo_name: mpleinti/pyspi

reason: import numpy,import astropy,from astropy


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

import numpy as np

from ard_gsm.qchem import QChem, QChemError
from ard_gsm.mol import MolGraph, SanitizationError
from ard_gsm.reaction import Reaction, group_reactions_by_products, group_reactions_by_connection_changes
from ard_gsm.util import write_xyz_file

def parse_reaction(reactant, prod_file, ts_file,
                   keep_isomorphic=False, edist_max=5.0, gdist_max=1.0,
                   normal_mode_check=False, soft_check=False, negative_barrier_check=True):
    """
    Parse reaction given the reactant MolGraph, the product log file,
    and the TS log file. Return None if there was an error in the TS
    job, the job was not valid, or other criteria are not met.
    """
    product = qchem2molgraph(prod_file, return_none_on_err=True, freq_only=True, print_msg=False)
    if product is None:
        print(f'Ignored {prod_file} because of negative frequency or error')
        return None

    if not keep_isomorphic and reactant.is_isomorphic(product):
        print(f'Ignored {prod_file} because product is isomorphic with reactant')
        return None

    ts_qts = qchem2molgraph(ts_file, return_qobj=True, return_none_on_err=True,
                            edist_max=edist_max, gdist_max=gdist_max, ts=True)
    if ts_qts is None:
        return None
    ts, qts = ts_qts

    rxn = Reaction(reactant, product, ts, product_file=prod_file, ts_file=ts_file)

    # Negative barriers shouldn't occur because we're calculating them
    # based on reactant/product wells, but check this just in case
    if negative_barrier_check:
        if ts.energy - reactant.energy < 0.0:
            print(f'Ignored {ts_file} because of negative barrier')
            return None
        elif ts.energy - product.energy < 0.0:
            print(f'Ignored {ts_file} because of negative reverse barrier')
            return None

    if normal_mode_check:
        normal_mode = qts.get_normal_modes()[0]  # First one corresponds to imaginary frequency
        if not rxn.normal_mode_analysis(normal_mode, soft_check=soft_check):
            print(f'Ignored {ts_file} because of failed normal mode analysis')
            return None

    return rxn

def remove_duplicates(reactions, ndup=1, group_by_connection_changes=False, atommap=True, set_smiles=True):
    """
    Group all reactions and remove all duplicates in each group keeping
    only at most ndup identical reactions. Also set the product SMILES.
    """
    extracted_reactions = {}

    # Group duplicate reactions
    if group_by_connection_changes:
        reaction_groups = group_reactions_by_connection_changes(reactions)
    else:
        reaction_groups = group_reactions_by_products(reactions)

    # Extract the ndup lowest barrier reactions from each group
    for group in reaction_groups:
        barriers = [(num, rxn.ts.energy - rxn.reactant.energy) for num, rxn in group.items()]
        barriers.sort(key=lambda x: x[1])

        nextracted = 0
        for num, _ in barriers:
            rxn = group[num]
            if set_smiles:
                try:
                    rxn.product_smiles = rxn.product.perceive_smiles(atommap=atommap)
                except SanitizationError:
                    print(f'Ignored number {num} because Smiles perception failed')
                    continue

            extracted_reactions[num] = rxn
            nextracted += 1

            if nextracted >= ndup:
                break

    return extracted_reactions

def qchem2molgraph(logfile, return_qobj=False, return_none_on_err=False, **kwargs):
    """
    Convert a Q-Chem logfile to a MolGraph object. Return the QChem
    object in addition to the MolGraph if return_qobj is True. Catch
    QChemError if return_none_on_err is True and return None. Options
    in kwargs are passed to valid_job. If the job is not valid, also
    return None.
    """
    try:
        q = QChem(logfile=logfile)
    except QChemError as e:
        if return_none_on_err:
            print(e)
            return None
        else:
            raise

    if not valid_job(q, **kwargs):
        return None

    energy = q.get_energy() + q.get_zpe()  # With ZPE
    symbols, coords = q.get_geometry()
    mol = MolGraph(symbols=symbols, coords=coords, energy=energy)
    mol.infer_connections()

    if return_qobj:
        return mol, q
    else:
        return mol

def valid_job(q, edist_max=None, gdist_max=None, ts=False, freq_only=False, print_msg=True):
    try:
        freqs = q.get_frequencies()
    except QChemError as e:
        print(e)
        return False
    nnegfreq = sum(1 for freq in freqs if freq < 0.0)

    if (ts and nnegfreq != 1) or (not ts and nnegfreq != 0):
        if print_msg:
            print(f'{q.logfile}: {nnegfreq} negative frequencies!')
        return False

    if freq_only:
        return True

    assert edist_max is not None and gdist_max is not None

    edist = abs(q.get_energy() - q.get_energy(first=True)) * 627.5095
    geo = q.get_geometry()[1]
    gdiff = geo.flatten() - q.get_geometry(first=True)[1].flatten()
    gdist = np.sqrt(np.dot(gdiff, gdiff) / len(geo))

    if edist > edist_max:
        if print_msg:
            print(f'{q.logfile}: Large energy change of {edist:.2f} kcal/mol!')
        return False
    if gdist > gdist_max:
        if print_msg:
            print(f'{q.logfile}: Large geometry change of {gdist:.2f} Angstrom!')
        return False

    return True

def rxn2xyzfile(rxn, path):
    symbols = (mol.get_symbols() for mol in (rxn.reactant, rxn.ts, rxn.product))
    coords = (mol.get_coords() for mol in (rxn.reactant, rxn.ts, rxn.product))
    comments = [rxn.reactant_smiles, '', rxn.product_smiles]
    write_xyz_file(path, symbols, coords, comments=comments)

reason: import numpy


text:
"""
Usage:
  lm.py [--vecs=VECS] [--labels=LABELS] [--text=TEXT]
        [--task=TASK]
        [--deprobe_dir=DEPROBE_DIR]
        [--display_examples=DISPLAY_EXAMPLES]
        [--device=DEVICE]
        [--n=N]
        [--wandb]

Options:
  -h --help                     show this help message and exit
  --vecs=VECS                   input vectors file
  --labels=LABELS               labels file
  --text=TEXT                   text file
  --task=TASK                   task type. between word_ind, sen_len, task [default: task]
  --deprobe_dir=DEPROBE_DIR     directory where the amnesic_probing files are located.
  --display_examples=DISPLAY_EXAMPLES       number of examples to display [default: 10]
  --device=DEVICE               cpu, cuda:0, cuda:1, ... [default: cpu]
  --n=N                         number of training examples [default: 100000]
  --wandb                       log using wandb

"""

import glob
import json
import os
from collections import Counter

import numpy as np
import pandas as pd
import wandb
from docopt import docopt
from sklearn.utils import shuffle

from amnesic_probing.tasks.data_preparation import get_appropriate_data
from amnesic_probing.tasks.lm import eval_lm_performance, eval_task_performance
from amnesic_probing.tasks.utils import data_projection, read_files, get_lm_vals

_, tokenizer, out_embed, bias = get_lm_vals('bert-base-uncased')

def log_wandb(arguments):
    task_name = arguments['--deprobe_dir'].split('models/lm/')[1]
    task_type = task_name.split('/')[0]
    layer = arguments['--vecs'].split('/')[-1].split('.')[0]

    labels = arguments['--labels'].split('.')[0]
    data_orig = labels.split('data/')[1].split('/')[0]
    print(labels)
    print(data_orig)
    dataset = data_orig.split('_output', 1)[0]
    masking = data_orig.rsplit('_', 1)[1]

    config = dict(
        property=task_type,
        encoder='bert-base-uncased',
        dataset=dataset,
        masking=masking,
        layer=layer
    )

    wandb.init(
        name=task_type + '_eval_iterations',
        project="amnesic_probing",
        tags=["lm", "eval", task_type],
        config=config,
    )

def load_deprobing_params(in_file):
    with open(in_file, 'r') as f:
        meta = json.load(f)
    return meta

if __name__ == '__main__':
    arguments = docopt(__doc__)

    deprobe_dir = arguments['--deprobe_dir']
    if not os.path.isdir(deprobe_dir):
        assert 'Deprobing directory does not exists...'

    use_wandb = arguments['--wandb']
    if use_wandb:
        log_wandb(arguments)

    vecs_train, labels_train, sentences_train = read_files(arguments['--vecs'],
                                                           arguments['--labels'],
                                                           arguments['--text'], ignore_special_tokens=True)
    vecs_dev, labels_dev, sentences_dev = read_files(arguments['--vecs'].replace('train', 'dev'),
                                                     arguments['--labels'].replace('train', 'dev'),
                                                     arguments['--text'].replace('train', 'dev'),
                                                     ignore_special_tokens=True)

    task = arguments['--task']

    (x_train, y_train, words_train), (x_dev, y_dev, words_dev) = get_appropriate_data(task, vecs_train, labels_train,
                                                                                sentences_train,
                                                                                vecs_dev, labels_dev,
                                                                                sentences_dev)

    pos2ind = {p: i for i, p in enumerate(sorted(set([item for sublist in labels_train for item in sublist])))}

    print('number of classes', len(pos2ind))
    print('most common class', Counter(y_dev).most_common(1)[0][1] / float(len(y_dev)))

    meta = load_deprobing_params(deprobe_dir + '/meta.json')
    n_coords = int(meta['removed_directions'])
    n_classes = len(pos2ind)

    if use_wandb:
        wandb.run.summary['n_classes'] = len(pos2ind)
        wandb.run.summary['majority'] = Counter(y_dev).most_common(1)[0][1] / float(len(y_dev))
        wandb.run.summary['removed_directions'] = n_coords

    x_train, y_train = shuffle(x_train, y_train, random_state=0, n_samples=min(len(y_train), int(arguments['--n'])))
    device = arguments['--device']

    results_dic = {}
    for proj_file in glob.glob(deprobe_dir + '/P_*.npy'):
        P = np.load(proj_file)

        print('evaluating performance')
        proj_iter = int(proj_file.split('/P_')[1].split('.npy')[0])

        removed_directions = int((proj_iter + 1) * n_classes)
        # in case of 2 classes, each inlp iteration we remove a single direction
        if n_classes == 2:
            removed_directions /= 2

        print(removed_directions)
        lm_results = eval_lm_performance(tokenizer, out_embed, bias, x_dev, words_dev, P,
                                         n_coords=int(removed_directions), device=device)

        print('removing property from inputs')

        x_train_no_label = data_projection(x_train, P)
        x_dev_no_label = data_projection(x_dev, P)

        task_results = eval_task_performance(x_train, y_train, x_dev, y_dev, x_train_no_label, x_dev_no_label, task)

        all_results = {**lm_results, **task_results}
        for k, v in all_results.items():
            print(k, v)

        results_dic[proj_iter] = all_results

    results_df = pd.DataFrame(results_dic).T
    results_df.index = results_df.index.set_names('iter')
    results_df.columns = results_df.columns
    results_df = results_df.reset_index()

    if use_wandb:
        table_data = results_df.values.tolist()
        wandb.log({"results": wandb.Table(data=table_data, columns=results_df.columns.tolist())})


text:
from dataset import PoST, val_dataset
from torch.utils.data import DataLoader
from models.model import PoSTModel

import tqdm
import os
import numpy as np

data_factory = {
    "PoST": PoST.PoST
}

def get_save_name(fname):
    save_name = os.path.basename(fname)
    save_name = ''.join(save_name.split('.')[:-1]) + '.txt'
    return save_name

def main(opt):
    if opt.data in data_factory:
        data = data_factory[opt.data]
        data = data(opt.data_root)

    val_data = val_dataset.ValDataset(data, opt)
    val_dataloader = DataLoader(val_data,
                                batch_size=1,
                                num_workers=opt.num_worker,
                                shuffle=False,
                                collate_fn=lambda x: x[0])
    model = PoSTModel(opt)

    if os.path.exists(opt.save_root):
        print(f"[WARNING] save root '{opt.save_root}' already exists. results will be replaced")

    for i, (imgs, annos, init, idx, others) \
                            in enumerate(val_dataloader):
        if 'seq' in others:
            seq = others['seq']
        else:
            raise ValueError('sequence name is not identified')

        img_files = others['img_files']

        outputs = {}
        save_dir = os.path.join(opt.save_root, seq)

        print(f'[VAL-{seq}] start inference')
        for j, img in enumerate(tqdm.tqdm(imgs)):
            if j == 0:
                print(f'[VAL-{seq}] initialize')
                cnt = model.initialize(imgs[j], init)
            else:
                cnt = model.inference(img)

            if len(idx) > 0:
                cnt = cnt[idx]

            save_name = get_save_name(img_files[j])
            outputs[save_name] = cnt

        save_dir = os.path.join(opt.save_root, seq)
        os.makedirs(save_dir, exist_ok=True)

        print(f'[VAL-{seq}] save outputs')
        for save_name, cnt in outputs.items():
            save_path = os.path.join(save_dir, save_name)
            np.savetxt(save_path, cnt)

max_stars_repo_name: ghnam-ken/PoST

max_issues_repo_name: ghnam-ken/PoST

max_forks_repo_name: ghnam-ken/PoST


text:
import numpy as np
import pandas

from ..general_trace.datasource import GeneralTraceDataSource

# Example random-number trace pseudo-data-source
class RandomTraceDataSource(GeneralTraceDataSource):

    # Create the datasource
    def __init__(self):
        super().__init__('')

        # Make some random data
        data = np.zeros(shape=10000000, dtype=np.dtype([('trn_idx', 'u8'), ('cycles', 'u8'), ('stat', 'i8')]))

        for i in range(1, len(data)):
            s = data[i - 1][2] + np.random.randint(-3, 4)
            data[i] = (i, i * 4 + np.random.randint(0, 3), s)

        # Construct a dataframe
        self._pdf = pandas.DataFrame(data)

    # Overriding GeneralTraceDataSource to give GeneralTraceAdapter adapter access to raw data
    def pdf(self) -> pandas.DataFrame:
        return self._pdf


text:
import numpy as np
import os
import keras
from keras import regularizers, losses
from keras.models import Sequential, Model
from keras.layers import Input, Dense, Dropout, Reshape, BatchNormalization, Softmax, Concatenate
from keras.utils import plot_model
from vaemodel import *

sampleLen = 5
numUnits = 3
enc_denseLayerSizes = [20, 20]
enc_denseLayerActivations = ['relu', 'relu']
enc_dropouts = [0, 0]
enc_batchnorms = [False, False]
dec_denseLayerSizes = [20, 20]
dec_denseLayerActivations = ['relu', 'relu']
dec_dropouts = [0, 0]
dec_batchnorms = [False, False]
inf_layerSize = 20

myModel = multiVAE(sampleLen, numUnits, enc_denseLayerSizes, enc_denseLayerActivations, enc_dropouts, enc_batchnorms, dec_denseLayerSizes, dec_denseLayerActivations, dec_dropouts, dec_batchnorms, inf_layerSize)

myModel.createFullNetwork()
myModel.createGenModel()


text:
from unittest.mock import patch

import numpy as np

from lambeq import QuantumModel

@patch.multiple(QuantumModel, __abstractmethods__=set())
def test_normalise():
    model = QuantumModel()
    inputs = np.linspace(-10,10,21)
    normalised = model._normalise_vector(inputs)
    assert abs(normalised.sum()-1.0)<1e-8
    assert np.all(normalised>=0)

max_stars_repo_name: CQCL/lambeq

max_issues_repo_name: CQCL/lambeq

max_forks_repo_name: CQCL/lambeq


text:
import bpy
from PIL import ImageDraw
import os
import torch
import cv2
import math
import numpy as np
from math import radians
from PIL import Image, ImageFilter
from dataset.data_utils import process_viewpoint_label, TransLightning, resize_pad, random_crop
import torchvision.transforms as transforms
from model.resnet import resnet50
from model.vp_estimator import BaselineEstimator, Estimator

# Crop and Resize the image and paste it to the bg
def crop_resize_paste(im, bg, left, upper, right, lower):
    # crop the RGBA image according to alpha channel
    bbox = im.getbbox()
    im = im.crop(bbox)

    # resize & padding the rendering then paste on the bg
    w, h = im.size
    target_w, target_h = right - left, lower - upper
    ratio = min(float(target_w) / w, float(target_h / h))
    new_size = (int(w * ratio), int(h * ratio))
    im = im.resize(new_size, Image.BILINEAR)
    bg.paste(im, (left + (target_w - new_size[0]) // 2, upper + (target_h - new_size[1]) // 2))

# create a lamp with an appropriate energy
def makeLamp(lamp_name, rad):
    # Create new lamp data block
    lamp_data = bpy.data.lamps.new(name=lamp_name, type='POINT')
    lamp_data.energy = rad
    # modify the distance when the object is not normalized
    # lamp_data.distance = rad * 2.5

    # Create new object with our lamp data block
    lamp_object = bpy.data.objects.new(name=lamp_name, object_data=lamp_data)

    # Link lamp object to the scene so it'll appear in this scene
    scene = bpy.context.scene
    scene.objects.link(lamp_object)
    return lamp_object

def parent_obj_to_camera(b_camera):
    # set the parenting to the origin
    origin = (0, 0, 0)
    b_empty = bpy.data.objects.new("Empty", None)
    b_empty.location = origin
    b_camera.parent = b_empty

    scn = bpy.context.scene
    scn.objects.link(b_empty)
    scn.objects.active = b_empty
    return b_empty

def clean_obj_lamp_and_mesh(context):
    scene = context.scene
    objs = bpy.data.objects
    meshes = bpy.data.meshes
    for obj in objs:
        if obj.type == "MESH" or obj.type == 'LAMP':
            scene.objects.unlink(obj)
            objs.remove(obj)
    for mesh in meshes:
        meshes.remove(mesh)

def angles_to_matrix(angles):
    """Compute the rotation matrix from euler angles for a mini-batch"""
    azi = angles[:, 0]
    ele = angles[:, 1]
    rol = angles[:, 2]
    element1 = (torch.cos(rol) * torch.cos(azi) - torch.sin(rol) * torch.cos(ele) * torch.sin(azi)).unsqueeze(1)
    element2 = (torch.sin(rol) * torch.cos(azi) + torch.cos(rol) * torch.cos(ele) * torch.sin(azi)).unsqueeze(1)
    element3 = (torch.sin(ele) * torch.sin(azi)).unsqueeze(1)
    element4 = (-torch.cos(rol) * torch.sin(azi) - torch.sin(rol) * torch.cos(ele) * torch.cos(azi)).unsqueeze(1)
    element5 = (-torch.sin(rol) * torch.sin(azi) + torch.cos(rol) * torch.cos(ele) * torch.cos(azi)).unsqueeze(1)
    element6 = (torch.sin(ele) * torch.cos(azi)).unsqueeze(1)
    element7 = (torch.sin(rol) * torch.sin(ele)).unsqueeze(1)
    element8 = (-torch.cos(rol) * torch.sin(ele)).unsqueeze(1)
    element9 = (torch.cos(ele)).unsqueeze(1)
    return torch.cat((element1, element2, element3, element4, element5, element6, element7, element8, element9), dim=1)

def render_obj(obj, output_dir, azi, ele, rol, name, shape=[512, 512], forward=None, up=None):
    clean_obj_lamp_and_mesh(bpy.context)

    # Setting up the environment
    scene = bpy.context.scene
    scene.render.resolution_x = shape[1]
    scene.render.resolution_y = shape[0]
    scene.render.resolution_percentage = 100
    scene.render.alpha_mode = 'TRANSPARENT'

    # Camera setting
    cam = scene.objects['Camera']
    cam_constraint = cam.constraints.new(type='TRACK_TO')
    cam_constraint.track_axis = 'TRACK_NEGATIVE_Z'
    cam_constraint.up_axis = 'UP_Y'
    b_empty = parent_obj_to_camera(cam)
    cam_constraint.target = b_empty

    # Light setting
    lamp_object = makeLamp('Lamp1', 5)
    lamp_add = makeLamp('Lamp2', 1)

    if forward is not None and up is not None:
        bpy.ops.import_scene.obj(filepath=obj, axis_forward=forward, axis_up=up)
    else:
        bpy.ops.import_scene.obj(filepath=obj)

    # normalize it and set the center
    for object in bpy.context.scene.objects:
        if object.name in ['Camera', 'Lamp'] or object.type == 'EMPTY':
            continue
        bpy.context.scene.objects.active = object
        max_dim = max(object.dimensions)
        object.dimensions = object.dimensions / max_dim if max_dim != 0 else object.dimensions

    # Output setting
    scene.render.image_settings.file_format = 'PNG'
    scene.render.filepath = os.path.join(output_dir, name)

    # scene.render.filepath = os.path.join(output_dir, name + '_render_%03d_%03d_%03d' % (int(azi), int(ele), int(rol)))

    # transform Euler angles from degrees into radians
    azi = radians(azi)
    ele = radians(ele)
    rol = radians(rol)

    r = 3
    loc_y = r * math.cos(ele) * math.cos(azi)
    loc_x = r * math.cos(ele) * math.sin(azi)
    loc_z = r * math.sin(ele)
    print(loc_x, loc_y, loc_z)

    print("===========")
    cam.location = (loc_x, loc_y, loc_z + 0.5)
    bpy.ops.render.render(write_still=True)

    cam.location = (loc_x, loc_y, loc_z + 0.5)
    lamp_object.location = (loc_x, loc_y, 10)
    lamp_add.location = (loc_x, loc_y, -10)

    bpy.ops.render.render(write_still=True)

    # Change the in-plane rotation
    cam_ob = bpy.context.scene.camera
    bpy.context.scene.objects.active = cam_ob  # select the camera object
    distance = np.sqrt(loc_x ** 2 + loc_y ** 2 + loc_z ** 2)
    bpy.ops.transform.rotate(value=90, axis=(loc_x / distance, loc_y / distance, loc_z / distance),
                             constraint_axis=(False, False, False), constraint_orientation='GLOBAL', mirror=False,
                             proportional='DISABLED', proportional_edit_falloff='SMOOTH',
                             proportional_size=1)
    print(loc_x / distance, loc_y / distance, loc_z / distance)
    print(rol)
    print(rol * 180 / math.pi)
    bpy.ops.render.render(write_still=True)

img_name = "../test_img/Screenshot 2022-01-26 120213.png"
img_name = "../test_img/horse.png"

obj_path = "../data/Pix3D/model/chair/IKEA_BERNHARD/model.obj"
save_name = "PoseImg.png"
out_dir = ".../test_img/out"

img = Image.open(img_name).convert('RGB')
im_pos = img.copy()

normalize = transforms.Normalize(mean=[0.485, 0.456, 0.406], std=[0.229, 0.224, 0.225])
im_transform = transforms.Compose([transforms.ToTensor(), normalize])

#im = im.crop((left, upper, right, lower))
im_flip = img.transpose(Image.FLIP_LEFT_RIGHT)
im = resize_pad(img, 224)
im_flip = resize_pad(im_flip, 224)
im = im_transform(im)
im_flip = im_transform(im_flip)
#im = im_flip

# ================CREATE NETWORK============================ #
os.environ["CUDA_VISIBLE_DEVICES"] = "0"

net_feat = resnet50(num_classes=128)
net_vp = BaselineEstimator(img_feature_dim=2048)

net_feat.cuda()
net_vp.cuda()

state = torch.load("../exps/PoseContrast_ObjectNet3d_ZeroShot/ckpt.pth")
net_feat.load_state_dict(state['net_feat'])
net_vp.load_state_dict(state['net_vp'])
net_feat.eval()
net_vp.eval()

im = im[None, :].cuda()
feat, _ = net_feat(im)
out = net_vp(feat)
vp_pred, score = net_vp.compute_vp_pred(out, True)

def RotMat(axis, ang):
    s = np.sin(ang)
    c = np.cos(ang)
    res = np.zeros((3, 3))
    if axis == 'Z':
        res[0, 0] = c
        res[0, 1] = -s
        res[1, 0] = s
        res[1, 1] = c
        res[2, 2] = 1
    elif axis == 'Y':
        res[0, 0] = c
        res[0, 2] = s
        res[1, 1] = 1
        res[2, 0] = -s
        res[2, 2] = c
    elif axis == 'X':
        res[0, 0] = 1
        res[1, 1] = c
        res[1, 2] = -s
        res[2, 1] = s
        res[2, 2] = c
    return res

def angle2dcm(azi, ele, theta):
    azi = radians(azi)
    ele = radians(ele)
    theta = radians(theta)
    PI = np.arccos(-1)
    return np.dot(RotMat('Z', theta), np.dot(RotMat('X', - (PI / 2 - ele)), RotMat('Z', - azi)))

#vp_pred[:, 1] = vp_pred[:, 1] - 180.
#vp_pred[:, 2] = vp_pred[:, 2] - 180.
#vp_pred = vp_pred * np.pi / 180.

    # render viewpoint predictions
[azi, ele, inp] = vp_pred[0].detach().cpu().numpy()

mat = angle2dcm(azi, ele - 90, inp - 180)
mat = torch.FloatTensor(mat)
front = torch.FloatTensor([0,0,1])
right = torch.FloatTensor([0,1,0])

up = torch.FloatTensor([1,0,0])
print(mat)
print("....")
print("Front", torch.matmul(front, mat))
print("Up", torch.matmul(up, mat))
print("Right", torch.matmul(right, mat))

render_obj(obj_path,
           "out",
           azi, ele - 90, inp - 180,
           "test.png",
           [521, 512], None, None)

#render_obj(obj_path,
#            "out",
#            azi, ele, inp,
#            "test.png",
#            [521, 512], None, None)

# blend viewpoint results into the original image
#img = cv2.imread(img_name)
#w, h = img.shape[1], img.shape[0]
#bg = Image.new("RGBA", (w, h))
#vp_rend = Image.open(os.path.join(out_dir, save_name))
#crop_resize_paste(vp_rend, bg, left, upper, right, lower)

#draw = ImageDraw.Draw(bg)
#draw.rectangle([(left, upper), (right, lower)], outline=(0, 255, 0, 255),
#               width=int(0.02 * min(img.shape[0], img.shape[1])))
#print("bgsave")
#bg.save(os.path.join(out_dir, save_name))

max_stars_repo_name: henlein/PoseContrast

max_issues_repo_name: henlein/PoseContrast

max_forks_repo_name: henlein/PoseContrast


text:
import collections
import copy
import os
import platform
from typing import Dict, Union, List, Any, Optional, Tuple

import numpy as np
import gym
import gym.spaces
from cffi import FFI

SCRIPT_DIR = os.path.dirname(os.path.abspath(__file__))
ACTION_KEY = "action"
STATE_NEEDS_RESET = "needs_reset"
STATE_WAIT_ACT = "wait_act"
STATE_WAIT_WAIT = "step_wait"

class CVecEnv:
    """
    An environment instance created by an EnvLib, uses the VecEnv interface.

    https://github.com/openai/baselines/blob/master/baselines/common/vec_env/__init__.py

    Args:
        num_envs: number of environments to create
        lib_dir: a folder containing either lib{name}.so (Linux), lib{name}.dylib (Mac), or {name}.dll (Windows)
        lib_name: name of the library (minus the lib part)
        c_func_defs: list of cdefs that are passed to FFI in order to define custom functions that can then be called with env.call_func()
        options: options to pass to the libenv_make() call for this environment
        debug: if set to True, check array data to make sure it matches the provided spaces
        reuse_arrays: reduce allocations by using the same numpy arrays for each reset(), step(), and render() call
    """

    def __init__(
        self,
        num_envs: int,
        lib_dir: str,
        lib_name: str = "env",
        c_func_defs: Optional[List[str]] = None,
        options: Optional[Dict] = None,
        debug: bool = False,
        reuse_arrays: bool = False,
    ) -> None:
        self._debug = debug
        self._reuse_arrays = reuse_arrays

        if options is None:
            options = {}
        options = copy.deepcopy(options)

        if platform.system() == "Linux":
            lib_filename = f"lib{lib_name}.so"
        elif platform.system() == "Darwin":
            lib_filename = f"lib{lib_name}.dylib"
        elif platform.system() == "Windows":
            lib_filename = f"{lib_name}.dll"
        else:
            raise Exception(f"unrecognized platform {platform.system()}")

        if c_func_defs is None:
            c_func_defs = []

        # load cdef for libenv.h
        libenv_cdef = ""
        with open(os.path.join(SCRIPT_DIR, "libenv.h")) as f:
            inside_cdef = False
            for line in f:
                if line.startswith("// BEGIN_CDEF"):
                    inside_cdef = True
                elif line.startswith("// END_CDEF"):
                    inside_cdef = False
                elif line.startswith("#if") or line.startswith("#endif"):
                    continue

                if inside_cdef:
                    line = line.replace("LIBENV_API", "")
                    libenv_cdef += line

        self._ffi = FFI()
        self._ffi.cdef(libenv_cdef)

        for cdef in c_func_defs:
            self._ffi.cdef(cdef)

        self._lib_path = os.path.join(lib_dir, lib_filename)
        assert os.path.exists(self._lib_path), f"lib not found at {self._lib_path}"
        # unclear if this is necessary, but nice to not have symbols conflict if possible
        dlopen_flags = (
            self._ffi.RTLD_NOW | self._ffi.RTLD_LOCAL  # pylint: disable=no-member
        )
        if platform.system() == "Linux":
            dlopen_flags |= self._ffi.RTLD_DEEPBIND  # pylint: disable=no-member
        self._c_lib = self._ffi.dlopen(name=self._lib_path, flags=dlopen_flags)
        # dlclose will be called automatically when the library goes out of scope
        # https://cffi.readthedocs.io/en/latest/cdef.html#ffi-dlopen-loading-libraries-in-abi-mode
        # on mac os x, the library may not always be unloaded when you expect
        # https://developer.apple.com/videos/play/wwdc2017/413/?time=1776
        # loading/unloading the library all the time can be slow
        # it may be useful to keep a reference to an environment (and thus the c_lib object)
        # to avoid this happening

        self._options = options
        self._state = STATE_NEEDS_RESET

        c_options, self._options_keepalives = self._convert_options(
            self._ffi, self._c_lib, options
        )

        self._c_env = self._c_lib.libenv_make(num_envs, c_options[0])

        self.reward_range = (float("-inf"), float("inf"))
        self.spec = None
        self.num_envs = num_envs
        self.observation_space = self._get_spaces(self._c_lib.LIBENV_SPACES_OBSERVATION)
        self._action_space = self._get_spaces(self._c_lib.LIBENV_SPACES_ACTION)
        self._info_space = self._get_spaces(self._c_lib.LIBENV_SPACES_INFO)
        self._render_space = self._get_spaces(self._c_lib.LIBENV_SPACES_RENDER)

        # allocate buffers
        self._observations, self._observation_buffers = self._allocate_dict_space(
            self.num_envs, self.observation_space
        )

        # we only use dict spaces for consistency, but action is always a single space
        # the private version is the dict space, while the public version is a single space
        assert len(self._action_space.spaces) == 1, "action space can only be 1 element"
        assert list(self._action_space.spaces.keys())[0] == ACTION_KEY
        self.action_space = self._action_space.spaces[ACTION_KEY]
        dict_actions, self._action_buffers = self._allocate_dict_space(
            self.num_envs, self._action_space
        )
        self._actions = dict_actions[ACTION_KEY]

        self._renders, self._renders_buffers = self._allocate_dict_space(
            self.num_envs, self._render_space
        )

        self.metadata = {"render.modes": list(self._render_space.spaces.keys())}

        self._infos, self._infos_buffers = self._allocate_dict_space(
            self.num_envs, self._info_space
        )

        self._rews, self._rews_buffer = self._allocate_array(
            self.num_envs, np.dtype("float32")
        )
        self._dones, self._dones_buffer = self._allocate_array(
            self.num_envs, np.dtype("bool")
        )
        assert np.dtype("bool").itemsize == 1

        c_step = self._ffi.new("struct libenv_step *")
        c_step.obs = self._observation_buffers
        # cast the pointer to the buffer to avoid a warning from cffi
        c_step.rews = self._ffi.cast(
            self._ffi.typeof(c_step.rews).cname, self._rews_buffer
        )
        c_step.dones = self._ffi.cast(
            self._ffi.typeof(c_step.dones).cname, self._dones_buffer
        )
        c_step.infos = self._infos_buffers
        self._c_step = c_step

        self.closed = False
        self.viewer = None

    def __repr__(self):
        return f"<CVecEnv lib_path={self._lib_path} options={self._options}>"

    def _numpy_aligned(self, shape, dtype, align=64):
        """
        Allocate an aligned numpy array, based on https://github.com/numpy/numpy/issues/5312#issuecomment-299533915
        """
        n_bytes = np.prod(shape) * dtype.itemsize
        arr = np.zeros(n_bytes + (align - 1), dtype=np.uint8)
        data_align = arr.ctypes.data % align
        offset = 0 if data_align == 0 else (align - data_align)
        view = arr[offset : offset + n_bytes].view(dtype)
        return view.reshape(shape)

    def _allocate_dict_space(
        self, num_envs: int, dict_space: gym.spaces.Dict
    ) -> Tuple[collections.OrderedDict, Any]:
        """
        Allocate arrays for a space, returns an OrderedDict of numpy arrays along with a backing bytearray
        """
        result = collections.OrderedDict()  # type: collections.OrderedDict
        length = len(dict_space.spaces) * num_envs
        buffers = self._ffi.new(f"void *[{length}]")
        for space_idx, (name, space) in enumerate(dict_space.spaces.items()):
            actual_shape = (num_envs,) + space.shape
            arr = self._numpy_aligned(shape=actual_shape, dtype=space.dtype)
            result[name] = arr
            for env_idx in range(num_envs):
                buffers[space_idx * num_envs + env_idx] = self._ffi.from_buffer(
                    arr.data[env_idx:]
                )
        return result, buffers

    def _allocate_array(self, num_envs: int, dtype: np.dtype) -> Tuple[np.ndarray, Any]:
        arr = self._numpy_aligned(shape=(num_envs,), dtype=dtype)
        return arr, self._ffi.from_buffer(arr.data)

    @staticmethod
    def _convert_options(ffi: Any, c_lib: Any, options: Dict) -> Any:
        """
        Convert a dictionary to libenv_options
        """
        keepalives = (
            []
        )  # add variables to here to keep them alive after this function returns
        c_options = ffi.new("struct libenv_options *")
        c_option_array = ffi.new("struct libenv_option[%d]" % len(options))

        for i, (k, v) in enumerate(options.items()):
            name = str(k).encode("utf8")
            assert (
                len(name) < c_lib.LIBENV_MAX_NAME_LEN - 1
            ), "length of options key is too long"
            if isinstance(v, bytes):
                c_data = ffi.new("char[]", v)
                dtype = c_lib.LIBENV_DTYPE_UINT8
                count = len(v)
            elif isinstance(v, str):
                c_data = ffi.new("char[]", v.encode("utf8"))
                dtype = c_lib.LIBENV_DTYPE_UINT8
                count = len(v)
            elif isinstance(v, bool):
                c_data = ffi.new("uint8_t*", v)
                dtype = c_lib.LIBENV_DTYPE_UINT8
                count = 1
            elif isinstance(v, int):
                assert -2 ** 31 < v < 2 ** 31
                c_data = ffi.new("int32_t*", v)
                dtype = c_lib.LIBENV_DTYPE_INT32
                count = 1
            elif isinstance(v, float):
                c_data = ffi.new("float*", v)
                dtype = c_lib.LIBENV_DTYPE_FLOAT32
                count = 1
            elif isinstance(v, np.ndarray):
                c_data = ffi.new("char[]", v.tobytes())
                if v.dtype == np.dtype("uint8"):
                    dtype = c_lib.LIBENV_DTYPE_UINT8
                elif v.dtype == np.dtype("int32"):
                    dtype = c_lib.LIBENV_DTYPE_INT32
                elif v.dtype == np.dtype("float32"):
                    dtype = c_lib.LIBENV_DTYPE_FLOAT32
                else:
                    assert False, f"unsupported type {v.dtype}"
                count = v.size
            else:
                assert False, f"unsupported value {v} for option {k}"

            c_option_array[i].name = name
            c_option_array[i].dtype = dtype
            c_option_array[i].count = count
            c_option_array[i].data = c_data
            keepalives.append(c_data)

        keepalives.append(c_option_array)
        c_options.items = c_option_array
        c_options.count = len(options)
        return c_options, keepalives

    @staticmethod
    def _check_arrays(
        arrays: Dict[str, np.ndarray], num_envs: int, dict_space: gym.spaces.Dict
    ) -> None:
        """
        Make sure each array in arrays matches the given dict_space
        """
        for name, space in dict_space.spaces.items():
            arr = arrays[name]
            assert isinstance(arr, np.ndarray)
            expected_shape = (num_envs,) + space.shape
            assert (
                arr.shape == expected_shape
            ), f"array is invalid shape expected={expected_shape} received={arr.shape}"
            assert (
                arr.dtype == space.dtype
            ), f"array has invalid dtype expected={space.dtype} received={arr.dtype}"
            if isinstance(space, gym.spaces.Box):
                # we only support single low/high values
                assert (
                    len(np.unique(space.low)) == 1 and len(np.unique(space.high)) == 1
                )
                low = np.min(space.low)
                high = np.max(space.high)
            elif isinstance(space, gym.spaces.Discrete):
                low = 0
                high = space.n - 1
            else:
                assert False, "unrecognized space type"
            assert np.min(arr) >= low, (
                '"%s" has values below space lower bound, %f < %f'
                % (name, np.min(arr), low)
            )
            assert np.max(arr) <= high, (
                '"%s" has values above space upper bound, %f > %f'
                % (name, np.max(arr), high)
            )

    def _maybe_copy_ndarray(self, obj: np.ndarray) -> np.ndarray:
        """
        Copy a single numpy array if reuse_arrays is False,
        otherwise just return the object
        """
        if self._reuse_arrays:
            return obj
        else:
            return obj.copy()

    def _maybe_copy_dict(self, obj: Dict[str, Any]) -> Dict[str, Any]:
        """
        Copy a list of dicts of numpy arrays if reuse_arrays is False,
        otherwise just return the object
        """
        if self._reuse_arrays:
            return obj
        else:
            result = {}
            for name, arr in obj.items():
                result[name] = arr.copy()
            return result

    def _get_spaces(self, c_name: Any) -> gym.spaces.Dict:
        """
        Get a c space and convert to a gym space
        """
        count = self._c_lib.libenv_get_spaces(self._c_env, c_name, self._ffi.NULL)
        if count == 0:
            return gym.spaces.Dict([])

        c_spaces = self._ffi.new("struct libenv_space[%d]" % count)
        self._c_lib.libenv_get_spaces(self._c_env, c_name, c_spaces)

        # convert to gym spaces
        spaces = []
        for i in range(count):
            c_space = c_spaces[i]

            name = self._ffi.string(c_space.name).decode("utf8")
            shape = []
            for j in range(c_space.ndim):
                shape.append(c_space.shape[j])

            if c_space.dtype == self._c_lib.LIBENV_DTYPE_UINT8:
                dtype = np.dtype("uint8")
                low = c_space.low.uint8
                high = c_space.high.uint8
            elif c_space.dtype == self._c_lib.LIBENV_DTYPE_INT32:
                dtype = np.dtype("int32")
                low = c_space.low.int32
                high = c_space.high.int32
            elif c_space.dtype == self._c_lib.LIBENV_DTYPE_FLOAT32:
                dtype = np.dtype("float32")
                low = c_space.low.float32
                high = c_space.high.float32
            else:
                assert False, "unknown dtype"

            if c_space.type == self._c_lib.LIBENV_SPACE_TYPE_BOX:
                space = gym.spaces.Box(shape=shape, low=low, high=high, dtype=dtype)
            elif c_space.type == self._c_lib.LIBENV_SPACE_TYPE_DISCRETE:
                assert shape == [1], "discrete space must have scalar shape"
                assert low == 0 and high > 0, "discrete low/high bounds are incorrect"
                space = gym.spaces.Discrete(n=high + 1)
                space.dtype = dtype
            else:
                assert False, "unknown space type"
            spaces.append((name, space))
        # c spaces are aways a single-layer Dict space
        return gym.spaces.Dict(spaces)

    def reset(self) -> Dict[str, np.ndarray]:
        """
        Reset the environment and return the first observation
        """
        self._state = STATE_WAIT_ACT

        self._c_lib.libenv_reset(self._c_env, self._c_step)
        return self._maybe_copy_dict(self._observations)

    def step_async(self, actions: np.ndarray) -> None:
        """
        Asynchronously take an action in the environment, doesn't return anything.
        """
        assert self._state == STATE_WAIT_ACT
        self._state = STATE_WAIT_WAIT

        if self._debug:
            self._check_arrays({"action": actions}, self.num_envs, self._action_space)

        self._actions[:] = actions
        self._c_lib.libenv_step_async(self._c_env, self._action_buffers, self._c_step)

    def step_wait(
        self
    ) -> Tuple[Dict[str, np.ndarray], np.ndarray, np.ndarray, List[Dict[str, Any]]]:
        """
        Step the environment, returns (obs, rews, dones, infos)
        """
        assert self._state == STATE_WAIT_WAIT
        self._state = STATE_WAIT_ACT

        self._c_lib.libenv_step_wait(self._c_env)

        if self._debug:
            self._check_arrays(
                self._observations, self.num_envs, self.observation_space
            )

        infos = [{} for _ in range(self.num_envs)]  # type: List[Dict]
        for key, values in self._infos.items():
            for env_idx in range(self.num_envs):
                v = values[env_idx]
                if v.shape == (1,):
                    # extract scalar values
                    v = v[0]
                infos[env_idx][key] = v

        return (
            self._maybe_copy_dict(self._observations),
            self._maybe_copy_ndarray(self._rews),
            self._maybe_copy_ndarray(self._dones),
            infos,
        )

    def step(
        self, actions: np.ndarray
    ) -> Tuple[Dict[str, np.ndarray], np.ndarray, np.ndarray, List[Dict[str, Any]]]:
        """
        Step the environment, combines step_async() and step_wait() and makes this interface mostly compatible with the normal gym interface
        """
        self.step_async(actions)
        return self.step_wait()

    def render(self, mode: str = "human") -> Union[bool, np.ndarray]:
        """
        Render the environment.

        mode='human' tells the environment to render in some human visible way and
        returns True if the user user visible window created by the environment is still open

        Other modes are up to the environment, but end up returning a numpy array of some kind
        that will be tiled as if it were an image by this code.  Call get_images() instead if
        that is not the behavior you want.
        """
        if (
            mode == "human"
            and "human" not in self.metadata["render.modes"]
            and "rgb_array" in self.metadata["render.modes"]
        ):
            # fallback human mode
            viewer = self.get_viewer()
            self._render(mode="rgb_array")
            images = self._maybe_copy_ndarray(self._renders["rgb_array"])
            viewer.imshow(self._tile_images(images))
            return viewer.isopen

        isopen = self._render(mode=mode)
        if mode == "human":
            return isopen
        else:
            images = self._maybe_copy_ndarray(self._renders[mode])
            return self._tile_images(images)

    def _render(self, mode) -> Union[bool, np.ndarray]:
        """Internal render method, returns is_open and updates the buffers in self._render"""
        assert self._state == STATE_WAIT_ACT
        assert mode in self.metadata["render.modes"], "unsupported render mode"

        c_mode = self._ffi.new("char[]", mode.encode("utf8"))
        # only pass the render buffers for the selected mode
        space_idx = list(self._render_space.spaces.keys()).index(mode)
        render_buffers = self._renders_buffers[
            space_idx * self.num_envs : (space_idx + 1) * self.num_envs
        ]
        return self._c_lib.libenv_render(self._c_env, c_mode, render_buffers)

    def _tile_images(self, images: np.ndarray) -> np.ndarray:
        """Tile a set of NHWC images"""
        num_images, height, width, chans = images.shape
        width_images = int(np.ceil(np.sqrt(num_images)))
        height_images = int(np.ceil(float(num_images) / width_images))
        result = np.zeros(
            (height_images * height, width_images * width, chans), dtype=np.uint8
        )
        for col in range(width_images):
            for row in range(height_images):
                idx = row * width_images + col
                if idx >= len(images):
                    continue
                result[
                    row * height : (row + 1) * height,
                    col * width : (col + 1) * width,
                    :,
                ] = images[idx]
        return result

    def get_images(self) -> np.ndarray:
        """
        Get rendered images from the environments, if supported.

        The returned array's shape will be (num_envs, height, width, num_colors)
        """
        self._render(mode="rgb_array")
        return self._maybe_copy_ndarray(self._renders["rgb_array"])

    def get_viewer(self):
        """Get the viewer instance being used by render()"""
        if self.viewer is None:
            from gym.envs.classic_control import rendering

            self.viewer = rendering.SimpleImageViewer()
        return self.viewer

    def close_extras(self):
        """
        Override this to close environment-specific resources without having to override close()'

        This method is guaranteed to only be called once, unlike close()
        """

    def close(self) -> None:
        """Close the environment and free any resources associated with it"""
        if not hasattr(self, "closed") or self.closed:
            return

        self.closed = True
        self._c_lib.libenv_close(self._c_env)
        self._c_lib = None
        self._ffi = None
        self._options_keepalives = None
        if self.viewer is not None:
            self.viewer.close()
        self.close_extras()

    def call_func(self, name: str, *args: Any) -> Any:
        """
        Call a function of the libenv declared in c_func_defs
        """
        return getattr(self._c_lib, name)(*args)

    @property
    def unwrapped(self):
        if hasattr(self, "venv"):
            return self.venv.unwrapped  # pylint: disable=no-member
        else:
            return self

    def seed(self, seed=None, env_idx=-1):
        """
        Seed the environment, this isn't used by VecEnvs but is part of the Gym Env API
        """
        self.reset_start_level(seed, env_idx)

    def reset_start_level(self, level_seed, env_idx=-1):
        self._c_lib.libenv_reset_start_level(self._c_env, level_seed, env_idx)

    def reset_at_index(self, env_idx) -> Dict[str, np.ndarray]:
        """
        Reset the environment and return the first observation
        """
        self._state = STATE_WAIT_ACT

        self._c_lib.libenv_reset_at_index(self._c_env, self._c_step, env_idx)
        return self._maybe_copy_dict(self._observations)

    def observe(self) -> Dict[str, np.ndarray]:
        """
        Reset the environment and return the first observation
        """
        self._state = STATE_WAIT_ACT

        self._c_lib.libenv_observe(self._c_env, self._c_step)
        return self._maybe_copy_dict(self._observations)

    def __del__(self):
        self.close()

max_stars_repo_name: minqi/procgen

max_issues_repo_name: minqi/procgen

max_forks_repo_name: minqi/procgen


text:
## buildMission20NC.py - Build Mission 20 Cleaned netCDF
# Import statements
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import netCDF4 as nc
from scipy import signal
from scipy import fft
from scipy import interpolate
from scipy import integrate
import cftime
import sys
import datetime
from halo import Halo

# Import DUNEX Tools
sys.path.append('..')
from tools import buildBasicMissionNC
from tools import missionView
from tools import microSWIFTTools

def main():
    '''
    @edwinrainville

    Description: Build and clean the mission 20 dataset.

    '''
    # Start loading spinnerfor first section
    spinner = Halo(text='Loading Dataset Information', spinner='dots')
    spinner.start()

    # Section 0:
    mission_num = 20

    # Define Project Directory
    project_dir = '../'

    # Define Data Directory
    data_dir = 'microSWIFT_data/'

    # Define Metadata Excel sheet name
    metadata_name = 'DUNEXMainExp_notes.xlsx'

    # Combine file name and project Directory
    metadata_filename = project_dir + metadata_name

    # Create dataframe object from DUNEX MetaData SpreadSheet
    dunex_xlsx = pd.read_excel(metadata_filename)

    # Get start and end times
    start_time = datetime.datetime.fromisoformat(dunex_xlsx['Start Time'].iloc[mission_num])
    end_time = datetime.datetime.fromisoformat(dunex_xlsx['End Time'].iloc[mission_num])

    # End Spinner for this section
    spinner.succeed('Dataset Information Loaded')
    spinner.stop()

    # Build a netCDF for the mission that has only raw time values read in
    spinner = Halo(text='Building Basic netCDF', spinner='dots')
    spinner.start()
    # mission_nc_path = buildBasicMissionNC.main(mission_num=mission_num)
    mission_nc_path = '../microSWIFT_data/cleanedDataset/mission_20.nc'
    spinner.succeed('Basic netCDF Built Successfully')
    spinner.stop()

    # Mask out invalid points in the dataset
    spinner = Halo(text='Masking Invalid Data Points in this Mission', spinner='dots')
    spinner.start()

    spinner.succeed('Invalid Data Masked Out Successfully')
    spinner.stop()

    # View the built mission netCDF
    spinner = Halo(text='Plotting Mission Data', spinner='dots')
    spinner.start()
    missionView.main(mission_nc_path, mission_num=mission_num)
    spinner.succeed('Mission Data Plotted Successfully')
    spinner.stop()

    # Compute

if __name__=='__main__':
    main()

max_stars_repo_name: edwinrainville/DUNEX

max_issues_repo_name: edwinrainville/DUNEXMainExp

max_forks_repo_name: edwinrainville/DUNEXMainExp

reason: import numpy,from scipy


text:
# -*- coding: utf-8 -*-
"""
Created on Mon May 17 18:53:35 2021

@author: rringuet
"""
from datetime import datetime, timezone
from numpy import vectorize
from os.path import isfile, basename
#read 1 day of data from cdf instead of from multiple .tec files

model_varnames={"Sigma_H":['Sigma_H','variable description',0,'SM','sph',['time','lon','lat'],"S"],
                "Sigma_P":['Sigma_P','variable description',1,'SM','sph',['time','lon','lat'],"S"],
                 "Phi_E":['Phi_E','variable description',2,'SM','sph',['time','lon','lat'],"W/m**2"],
                 "AveE_avgE":['E_avg','variable description',3,'SM','sph',['time','lon','lat'],'eV'],
                 "j_R":["j_R",'variable description',4,'SM','sph',['time','lon','lat'],"muA/m**2"],
                 "Phi":["Phi",'variable description',5,'SM','sph',['time','lon','lat'],"kV"],
                 "E_x":["E_x",'variable description',6,'SM','sph',['time','lon','lat'],"mV/m"],
                 "E_y":["E_y",'variable description',7,'SM','sph',['time','lon','lat'],"mV/m"],
                 "E_z":["E_z",'variable description',8,'SM','sph',['time','lon','lat'],"mV/m"],
                 "j_x":["j_x",'variable description',9,'SM','sph',['time','lon','lat'],"muA/m**2"],
                 "j_y":["j_y",'variable description',10,'SM','sph',['time','lon','lat'],"muA/m**2"],
                 "j_z":["j_z",'variable description',11,'SM','sph',['time','lon','lat'],"muA/m**2"],
                 "v_x":['v_x','variable description',12,'SM','sph',['time','lon','lat'],"km/s"],
                 "v_y":['v_y','variable description',13,'SM','sph',['time','lon','lat'],"km/s"],
                 "v_z":['v_z','variable description',14,'SM','sph',['time','lon','lat'],"km/s"],
                 "Q_Joule":['Q_Joule','variable description',15,'SM','sph',['time','lon','lat'],"mW/m**2"],
                 "Phi_nion":['Phi_nion','variable description',16,'SM','sph',['time','lon','lat'],"1/cm**2/s"],
                 "Binv_RT":['Binv_RT','variable description',17,'SM','sph',['time','lon','lat'],"1/T"],
                 "rho_RT":['rho_RT','variable description',18,'SM','sph',['time','lon','lat'],"amu/cm**3"],
                 "P_RT":['P_RT','variable description',19,'SM','sph',['time','lon','lat'],"Pa"],
                 "dLat_star":['dLat_star','variable description',20,'SM','sph',['time','lon','lat'],"deg"],
                 "dlon_star":['dlon_star','variable description',21,'SM','sph',['time','lon','lat'],"deg"]}

'''
#documentation variables:
"X":['x','km'],"Y":['y','km'],"Z":['z','km'],    (ignored, given in R_E on a unit sphere)
"Theta":['theta',"deg"],"Psi":['psi',"deg"],     (used as coordinates)
"Btilt_theta":['theta_Btilt',"deg"], "Btilt_psi":['psi_Btilt',"deg"]
(added directly to object for documentation purposes)
'''
@vectorize
def dts_to_hrs(datetime_string, filedate):
    '''Get hours since midnight from datetime string'''

    return (datetime.strptime(datetime_string, '%Y-%m-%d %H:%M:%S').replace(tzinfo=timezone.utc)\
            -filedate).total_seconds()/3600.

@vectorize
def filename_to_dts(filename, string_date):
    '''Get datetime string in format "YYYY-MM-SS HH:mm:SS" from filename'''

    mmhhss = basename(filename)[12:18]
    return string_date+' '+mmhhss[:2]+':'+mmhhss[2:4]+':'+mmhhss[4:]

def dts_to_ts(file_dts):
    '''Get datetime timestamp in UTC from datetime string'''

    return datetime.timestamp(datetime.strptime(file_dts, '%Y-%m-%d %H:%M:%S'
                                                ).replace(tzinfo=timezone.utc))

def ts_to_hrs(time_val, filedate):
    '''Convert utc timestamp to hours since midnight on filedate.'''

    return (datetime.utcfromtimestamp(time_val).replace(tzinfo=timezone.utc)-filedate).total_seconds()/3600.

#base for main class object----------------------------------------
#main class object
#file_dir = 'C:/Users/rringuet/Kamodo_WinDev1/SWMF_IE/'
#file_prefix = 'C:/Users/rringuet/Kamodo_WinDev1/SWMF_IE/i_e20180801'  #example
#files = glob.glob(file_dir+'i_e*.tec')  #for wrapper, this and the next line
#file_patterns = unique([file_dir+f.split('/')[-1].split('\\')[-1][:11] for f in files])
def MODEL():
    from numpy import array, NaN, abs, unique, append, zeros, diff, where, insert, flip
    from time import perf_counter
    from netCDF4 import Dataset
    from kamodo import Kamodo
    #print('KAMODO IMPORTED!')
    from kamodo_ccmc.readers.reader_utilities import regdef_3D_interpolators

    class MODEL(Kamodo):
        def __init__(self, full_file_prefix, variables_requested=[], runname="noname",
                     filetime=False, verbose=False, gridded_int=True, printfiles=False,
                     fulltime=True, **kwargs):
            '''file_prefix must be of form "3D***_tYYMMDD" to load all files for one day
             and include a complete path to the files'''
            super(MODEL, self).__init__()

            #check if given .nc file exists. If not, convert files with same prefix to netCDF
            file_prefix = basename(full_file_prefix)
            file_dir = full_file_prefix.split(file_prefix)[0]
            if not isfile(full_file_prefix+'.nc'):
                from kamodo_ccmc.readers.swmfie_tocdf import convertSWMFIE_toCDF
                test = convertSWMFIE_toCDF(full_file_prefix)
                if not test:
                    self.conversion_test = test  #only needed for 1 file/time cases
                    return    #if file conversion fails, return
                else: self.conversion_test = True
            else: self.conversion_test = True
            t0 = perf_counter()

            #determine type of prefix: for a day or for a hour
            if '-' in file_prefix: day_flag=False
            else: day_flag=True

            #establish time attributes first for file searching
            file_datestr = file_prefix[3:11]
            string_date = file_datestr[:4]+'-'+file_datestr[4:6]+'-'+file_datestr[6:8]  #'YYYY-MM-DD'
            self.filedate = datetime.strptime(string_date+' 00:00:00', \
                                              '%Y-%m-%d %H:%M:%S').replace(tzinfo=timezone.utc) #dt object

            #establish beginning and end time of file list
            cdf_data = Dataset(full_file_prefix+'.nc', 'r')
            files = cdf_data.file.split(',')
            self.datetimes = list(filename_to_dts([files[0], files[-1]], string_date))  #strings in format = YYYY-MM-DD HH:MM:SS
            self.filetimes=[dts_to_ts(file_dts) for file_dts in self.datetimes]   #timestamps in UTC
            t = array(cdf_data.variables['time'])  #hours since midnight
            if len(t)>1: self.dt = diff(t).max()*3600.  #t is in hours since midnight
            else: self.dt = 0

            if filetime and not fulltime: #(used when searching for neighboring files below)
                return  #return times as is to prevent recursion

            #if variables are given as integers, convert to standard names
            if len(variables_requested)>0:
                if isinstance(variables_requested[0], int):
                    tmp_var = [value[0] for key, value in model_varnames.items()\
                                           if value[2] in variables_requested]
                    variables_requested = tmp_var

            if fulltime:  #add boundary time (default value)
                #find other files with same pattern
                from glob import glob

                files = sorted(glob(file_dir+'i_e*'))
                if day_flag:
                    file_prefixes = unique([basename(f)[:11] for f in files\
                                            if '.nc' not in basename(f)])
                else:  #give prefix for hourly files
                    file_prefixes = unique([basename(f)[:14] for f in files\
                                            if '.nc' not in basename(f)])

                #find closest file by utc timestamp
                #swmf_ie has an open time at the end, so need a beginning time from the next file
                #files are automatically sorted by YYMMDD, so next file is next in the list
                current_idx = where(file_prefixes==file_prefix)[0]
                if current_idx+1==len(file_prefixes):
                    print('No later file available.')
                    filecheck = False
                    if filetime:
                        return
                else:
                    min_file_prefix = file_dir+file_prefixes[current_idx+1][0]  #+1 for adding an end time
                    kamodo_test = MODEL(min_file_prefix, filetime=True, fulltime=False)
                    if not kamodo_test.conversion_test:
                        print('No later file available.')
                        filecheck = False
                        if filetime:
                            return
                    else:
                        time_test = abs(kamodo_test.filetimes[0]-self.filetimes[1])
                        if time_test<=self.dt:  #if nearest file time at least within one timestep (hrs)
                            filecheck = True

                            #time only version if returning time for searching
                            if filetime:
                                kamodo_neighbor = MODEL(min_file_prefix, fulltime=False, filetime=True)
                                self.datetimes[1] = kamodo_neighbor.datetimes[0]
                                self.filetimes[1] = kamodo_neighbor.filetimes[0]
                                return  #return object with additional time (for SF code)

                            #get kamodo object with same requested variables to add to each array below
                            if verbose: print(f'Took {perf_counter()-t0:.3f}s to find closest file.')
                            kamodo_neighbor = MODEL(min_file_prefix,
                                                    variables_requested=variables_requested,
                                                    fulltime=False)
                            self.datetimes[1] = kamodo_neighbor.datetimes[0]
                            self.filetimes[1] = kamodo_neighbor.filetimes[0]
                            short_data = kamodo_neighbor.short_data
                            if verbose: print(f'Took {perf_counter()-t0:.3f}s to get data from closest file.')
                        else:
                            print(f'No later file found within {self.dt:.1f}s.')
                            filecheck = False
                            if filetime:
                                return

            #perform initial check on variables_requested list
            if len(variables_requested)>0 and fulltime:
                test_list = [value[0] for key, value in model_varnames.items()]
                err_list = [item for item in variables_requested if item not in test_list]
                if len(err_list)>0: print('Variable name(s) not recognized:', err_list)

            #get list of variables possible in these files using first file
            if len(variables_requested)>0:
                gvar_list = [key for key, value in model_varnames.items() \
                                 if value[0] in variables_requested and \
                                     key in cdf_data.variables.keys()]  # file variable names
                #check for variables requested but not available
                if len(gvar_list)!=len(variables_requested):
                    err_list = [value[0] for key, value in model_varnames.items() \
                                 if value[0] in variables_requested and \
                                     key not in cdf_data.variables.keys()]
                    if len(err_list)>0: print('Some requested variables are not available:', err_list)
            else:
                avoid_list = ['theta_Btilt', 'psi_Btilt']
                gvar_list = [key for key in cdf_data.variables.keys() \
                             if key in model_varnames.keys() and \
                                 key not in avoid_list]

            # Store variable's data and units, transposing the 2D+time array.
            variables = {model_varnames[key][0]:{'units':model_varnames[key][-1],
                                   'data':array(cdf_data.variables[key])}\
                              for key in gvar_list}
            self.theta_Btilt = array(cdf_data.variables['theta_Btilt'])
            self.psi_Btilt = array(cdf_data.variables['psi_Btilt'])

            #prepare and return data only for first timestamp
            if not fulltime:
                cdf_data.close()
                variables['time'] = self.filetimes[0]
                variables['theta_Btilt'] = self.theta_Btilt[0]
                variables['psi_Btilt'] = self.psi_Btilt[0]
                self.short_data = variables
                return

            #return if only one file found because interpolator code will break
            if len(files)<2:
                print('Not enough files found with given file prefix.')
                return

            #store variables
            self.filename = files
            self.runname = runname
            self.missing_value = NaN
            self.modelname = 'SWMF_IE'
            self._registered = 0
            if printfiles:
                print('Files:')
                for file in self.filename: print(file)

            #### Store coordinate data as class attributes
            if filecheck:
                new_time = ts_to_hrs(short_data['time'], self.filedate)  #new time in hours since midnight
                self._time = append(t, new_time)
            else:
                self._time = t

            #store coordinate data
            #self._radius = array(cdf_data.variables['radius'])
            self._lat = array(cdf_data.variables['lat'])
            self._lon = array(cdf_data.variables['lon'])
            cdf_data.close()
            if verbose: print(f'Took {perf_counter()-t0:.6f}s to read in data')

            #register interpolators for each variable
            varname_list, self.variables = [key for key in variables.keys()], {}  #store original list b/c gridded interpolators
            t_reg = perf_counter()
            for varname in varname_list:  #all are 3D variables
                if filecheck:  #if neighbor found
                    #append data for first time stamp, transpose and register
                    data_shape = list(variables[varname]['data'].shape)
                    data_shape[0]+=1  #add space for time
                    new_data = zeros(data_shape)
                    new_data[:-1,:,:] = variables[varname]['data']  #put in current data
                    new_data[-1,:,:] = short_data[varname]['data'][0,:,:]  #add in data for additional time
                else:
                    new_data = variables[varname]['data']
                self.variables[varname] = dict(units = variables[varname]['units'], data = new_data)
                self.register_3D_variable(self.variables[varname]['units'],
                                          self.variables[varname]['data'], varname,
                                          gridded_int)
            if verbose: print(f'Took {perf_counter()-t_reg:.5f}s to register '+\
                              f'{len(varname_list)} variables.')
            if verbose: print(f'Took a total of {perf_counter()-t0:.5f}s to kamodofy '+\
                              f'{len(gvar_list)} variables.')

        #define and register a 3D variable-----------------------------------------
        def register_3D_variable(self, units, variable, varname, gridded_int):
            """Registers a 3d interpolator with 3d signature"""

            #define and register the interpolators
            xvec_dependencies = {'time':'hr','lon':'deg','lat':'deg'}
            self = regdef_3D_interpolators(self, units, variable, self._time,
                                           self._lon, self._lat, varname,
                                           xvec_dependencies, gridded_int)
            return
    return MODEL

max_stars_repo_name: EnsembleGovServices/Kamodo-ccmc-readers

max_issues_repo_name: EnsembleGovServices/Kamodo-ccmc-readers

max_forks_repo_name: EnsembleGovServices/Kamodo-ccmc-readers

reason: from numpy


text:
# -*- coding: utf-8 -*-
from LLL import *
import numpy as np

# LLL
mat1 = LLLasMatCol()

v_1 = mat1[0]
v_2 = mat1[1]
v_3 = mat1[2]

print("Ejemplo : ")
print("A : ")
mat1 = create_matrix(mat1)
mat1 = lll_reduction(mat1)
print_mat(mat1)

print("¿A está reducida? ", islll(mat1))
print()

print(v_1)
print(v_2)
print(v_3)
print()

# GAUSS
mat2 = GaussMatCol()

u_1 = mat2[0]
u_2 = mat2[1]

print("Ejemplo : ")
print("A : ")
mat2 = create_matrix(mat2)
print_mat(mat2)
print("¿A está reducida? ", islll(mat2))
print()

print(u_1)
print(u_2)

max_stars_repo_name: JazzzFM/Latticed-Manim

max_issues_repo_name: JazzzFM/Latticed-Manim

max_forks_repo_name: JazzzFM/Latticed-Manim


text:
"""The :mod:`butterworth` module provides functions to design butterwort filters.

"""
from numpy import (mod, exp, zeros, ones, arange, kron, real, flipud,
                   conj, pi, fliplr, sqrt, tan, tile, concatenate, append, double)
from scipy.signal import butter, buttord, buttap

from pyfilterbank.sosfiltering import bilinear_sos

lowpass = 'lowpass'
highpass = 'highpass'
bandpass = 'bandpass'
bandstop = 'bandstop'

def butter_sos(band, L, v1, v2=0.0):
    """Compute weights of a digital Butterworth filter in cascade form.

    Parameters
    ----------
    band : {'lowpass', 'highpass', 'bandpass', 'bandstop'}
    L : int
        Number of lowpass poles. L is doubled for
        'bandpass' and 'bandstop'. L must be even.
    v1 : scalar
        First critical frequency (Hz-s); 0.0 <v1 < 0.5.
    v2 : scalar
        Second critical frequency; v1 < v2 < 0.5.
        v2 is used only if 'bandpass' or 'bandstop'.

    Returns
    -------
    sosmat : ndarray
        Contains the numerator and denominator coeffs for each
        cascade in one row.

    Notes
    -----
    Adapted from: Samuel D. Stearns, "Digital Signal Processing
    with Examples in MATLAB"
    """

    # Check for errors;
    if(v1<=0 or v1>=.5):
       raise Exception('Argument v1 must be >0.0 and <0.5')
    elif v2!=0 and (v2<=v1 or v2>=0.5):
       raise Exception('Argument v2 must be >v1 and <0.5')

   # Get the anlalog weights and convert to digital.
    if(v2 == 0):
       d, c = butter_analog_sos(band, L, tan(pi*v1))
    else:
       d, c = butter_analog_sos(band, L, tan(pi*v1), tan(pi*v2))
    b, a = bilinear_sos(d, c)
    # TODO: cut this step and give back b, a
    # but the other functions needing an sosmat should be adapted.
    sosmat = flipud(concatenate((b, a), axis=1))
    return sosmat

def butter_analog_sos(band, L, w1, w2=0):
    """Returns analog filter coeffitients for Butterworth filters.
    compute analog weights of a butterworth filter

    Parameters
    ----------
    band : {'lowpass', 'highpass, 'bandpass', 'bandstop'}
    L : int
        Order of lowpass / highpass filter.
    w1 : scalar
        Critical frequency one.
    w2 : scalar
        Critical frequency two. (for 'bandpass' or 'bandstop').

    Returns
    -------
    d, c :  Analog weights of the filter

    Notes
    -----
    implements SOS H(s)  butterwort
    if you need H(z) apply a bilinear transform

    """
    band = band.lower()
    L2 = int(L / 2.0)

    if mod(L, 2):
        raise Exception('Number of poles L must be even')
    if w1 <= 0:
        raise Exception('Frequency w1 must be in rad/s and >0')

    # define center frequency wc:
    if band==lowpass or band==highpass:
        wc = w1
    else:
        wc = w2-w1

    p = wc * exp(-1j * (2*arange(1,L2+1,dtype=double) + L-1) * pi / (2*L))

    # defining the lowpass filter:
    d = zeros((L2, 2), dtype=double).astype(complex)
    c = ones((L2, 2), dtype=double).astype(complex)
    d[:,1] = wc * ones(L2, dtype=double)
    c[:,1] = -p

    if band == highpass:
        d[:, 0] = d[:, 1]
        d[:, 1] = 0.0
        c = fliplr(c)
        c[:, 1] = c[:, 1] * wc**2

    elif band == bandpass:
        d[:, 0] = d[:, 1]
        d[:, 1] = zeros(L2)
        d = append(d, d, axis=0)
        d[L2:, 0] = zeros(L2)
        d[L2:, 1] = ones(L2)
        root = sqrt(c[:, 1]**2 - 4*c[:, 0]**2 * w1*w2)
        r1 = (-c[:, 1] + root) / (2 * c[:, 0])
        r2 = (-c[:, 1] - root) / (2 * c[:, 0])
        c[:, 0] = c[:, 0]
        c[:, 1] = -c[:, 0] * r1
        c = append(c, c, axis=0)
        c[L2:, 0] = 1.0
        c[L2:, 1] = -r2

    elif band == bandstop:
        root = sqrt(d[:, 0]**2 * wc**4 - 4*d[:, 1]**2 * w1*w2)
        r1 = (-d[:, 0] * wc**2 + root) / (2 * d[:, 1])
        r2 = (-d[:, 0] * wc**2 - root) / (2 * d[:, 1])
        d[:, 0] = d[:L2, 1]
        d[:, 1] = -d[:L2, 1] * r1
        d = append(d, d, axis=0)
        d[L2:2*L2, 0] = ones(L2)
        d[L2:2*L2, 1] = -r2
        root = sqrt(c[:, 0]**2 * wc**4 - 4*c[:, 1]**2 * w1*w2)
        r1 = (-c[:, 0] * wc**2 + root) / (2 * c[:, 1])
        r2 = (-c[:, 0] * wc**2 - root) / (2 * c[:, 1])
        c[:, 0] = c[:L2, 1]
        c[:, 1] = -c[:L2, 1] * r1
        c = append(c, c, axis=0)
        c[L2:2*L2, 0] = ones(L2)
        c[L2:2*L2, 1] = -r2

    return d, c

reason: from numpy,from scipy