text:
# -*- coding: utf-8 -*-
import numpy as np
import pyximport
import torch
from torch import nn

from bert.modeling import ACT2FN, BertLayerNorm, BertModel, BertPreTrainedModel
from trees import InternalParseNode, LeafParseNode

pyximport.install(
    build_dir=".caches", setup_args={"include_dirs": np.get_include()}, language_level=3
)

import chart_decoder

class ChartParser(BertPreTrainedModel):
    def __init__(
        self,
        config,
        tag_encoder,
        label_encoder,
        lstm_layers,
        lstm_dim,
        tag_embedding_dim,
        label_hidden_dim,
        dropout_prob,
    ):
        super(ChartParser, self).__init__(config)

        self.bert = BertModel(config)

        self.tag_encoder = tag_encoder
        self.label_encoder = label_encoder

        self.lstm = nn.LSTM(
            input_size=tag_embedding_dim + config.hidden_size,
            hidden_size=lstm_dim,
            num_layers=lstm_layers,
            batch_first=True,
            dropout=dropout_prob,
            bidirectional=True,
        )

        self.tag_embeddings = nn.Embedding(
            num_embeddings=len(tag_encoder),
            embedding_dim=tag_embedding_dim,
            padding_idx=tag_encoder.transform("[PAD]"),
        )

        self.hidden_dense = nn.Linear(
            in_features=3 * lstm_dim * 2, out_features=label_hidden_dim
        )

        self.label_classifier = nn.Linear(
            in_features=label_hidden_dim,
            out_features=len(label_encoder) - 1,  # Skip label ()
        )

        self.intermediate_act_fn = ACT2FN[config.hidden_act]

        self.layer_norm = BertLayerNorm(label_hidden_dim, eps=config.layer_norm_eps)

        self.dropout = nn.Dropout(dropout_prob)

        self.apply(self.init_bert_weights)

    def forward(self, ids, attention_masks, tags, sections, sentences, gold_trees=None):
        is_training = gold_trees is not None

        subtoken_embeddings, _ = self.bert(
            input_ids=ids,
            attention_mask=attention_masks,
            output_all_encoded_layers=False,
        )

        tag_embeddings = self.tag_embeddings(tags)

        sentence_sections = [0]
        span_embeddings = []

        # Loop over each sample in a mini-batch
        for _subtoken_embeddings, _tag_embeddings, section, sentence in zip(
            subtoken_embeddings, tag_embeddings, sections, sentences
        ):
            num_tokens = len(section)
            num_subtokens = sum(section)

            # Remove paddings
            _tag_embeddings = _tag_embeddings.narrow(dim=0, start=0, length=num_tokens)
            _subtoken_embeddings = _subtoken_embeddings.narrow(
                dim=0, start=1, length=num_subtokens
            )

            # Merge subtoken embeddings to form a single token embedding
            token_embeddings = []
            for _token_embeddings in _subtoken_embeddings.split(section, dim=0):
                token_embeddings.append(_token_embeddings.mean(dim=0, keepdim=True))

            token_embeddings = torch.cat(token_embeddings, dim=0)

            token_embeddings = torch.cat([_tag_embeddings, token_embeddings], dim=-1)

            token_embeddings = self.dropout(token_embeddings)

            lstm_embeddings, _ = self.lstm(token_embeddings.unsqueeze(dim=0))

            lstm_embeddings = lstm_embeddings.squeeze(dim=0)

            # Generate spans
            for length in range(1, len(sentence) + 1):
                for left in range(0, len(sentence) + 1 - length):
                    right = left + length

                    left_embedding = lstm_embeddings[left]

                    right_embedding = lstm_embeddings[right - 1]

                    average_embedding = lstm_embeddings.narrow(
                        dim=0, start=left, length=length
                    ).mean(dim=0)

                    span_embeddings.append(
                        torch.cat([left_embedding, average_embedding, right_embedding])
                    )

            sentence_sections.append(len(span_embeddings))

        span_embeddings = torch.stack(span_embeddings)

        label_scores = self.label_classifier(
            self.dropout(
                self.layer_norm(
                    self.intermediate_act_fn(self.hidden_dense(span_embeddings))
                )
            )
        )

        # Add score 0 for NULL label
        label_scores = torch.cat(
            (label_scores.new_zeros(label_scores.size(0), 1), label_scores), dim=-1
        )

        total_augmented_amount = 0.0

        all_pred_trees = []

        all_pred_indices = []
        all_pred_labels = []

        all_gold_indices = []
        all_gold_labels = []

        for sentence_idx, np_label_scores in enumerate(
            np.split(label_scores.detach().cpu().numpy(), sentence_sections[1:-1])
        ):
            decoder_args = {
                "sentence_len": len(sections[sentence_idx]),
                "num_previous_indices": sentence_sections[sentence_idx],
                "label_scores": np_label_scores,
                "is_training": is_training,
                "gold_tree": gold_trees and gold_trees[sentence_idx],
                "label_encoder": self.label_encoder,
            }

            _, pred_included_i, pred_included_j, pred_included_indices, pred_included_labels, pred_augmented_amount = chart_decoder.decode(
                False, **decoder_args
            )

            all_pred_indices.append(pred_included_indices)
            all_pred_labels.append(pred_included_labels)

            total_augmented_amount += pred_augmented_amount

            if is_training:
                _, _, _, gold_included_indices, gold_included_labels, _ = chart_decoder.decode(
                    True, **decoder_args
                )

                all_gold_indices.append(gold_included_indices)
                all_gold_labels.append(gold_included_labels)
            else:
                stack_idx = -1

                def make_tree():
                    nonlocal stack_idx
                    stack_idx += 1

                    pos_i, pos_j, label_index = (
                        pred_included_i[stack_idx],
                        pred_included_j[stack_idx],
                        pred_included_labels[stack_idx],
                    )

                    pos_i = pos_i.item()
                    pos_j = pos_j.item()

                    pred_label = self.label_encoder.inverse_transform(label_index)

                    if pos_i + 1 >= pos_j:
                        tag, word = sentences[sentence_idx][pos_i]
                        tree = LeafParseNode(pos_i, tag, word)
                        if pred_label:
                            tree = InternalParseNode(pred_label, [tree])
                        return [tree]
                    else:
                        left_trees = make_tree()
                        right_trees = make_tree()
                        children = left_trees + right_trees
                        if pred_label:
                            return [InternalParseNode(pred_label, children)]
                        else:
                            return children

                all_pred_trees.append(make_tree()[0])

        # Is training
        if is_training:
            all_pred_indices = ids.new_tensor(np.concatenate(all_pred_indices))
            all_pred_labels = ids.new_tensor(np.concatenate(all_pred_labels))

            all_gold_indices = ids.new_tensor(np.concatenate(all_gold_indices))
            all_gold_labels = ids.new_tensor(np.concatenate(all_gold_labels))

            loss = (
                label_scores[all_pred_indices, all_pred_labels].sum()
                - label_scores[all_gold_indices, all_gold_labels].sum()
                + total_augmented_amount
            )

            return loss

        return all_pred_trees

max_stars_repo_name: dnanhkhoa/bert-span-parser

max_issues_repo_name: dnanhkhoa/bert-span-parser

max_forks_repo_name: dnanhkhoa/bert-span-parser


text:
# -*- coding: utf-8 -*-

from keras import models
from keras import layers
from keras.datasets import mnist
from keras.utils import to_categorical
from keras import backend as K
from keras.preprocessing import image
import matplotlib.pyplot as plt
import numpy as np
import gc
import PIL

data_dir = '../DeepLearning-Learning_DATAS/03_FineTuning/'
h5_name = 'janken_small_cnn_01.h5'
#image_name = 'validation/gu/20180801_211326790991.jpg'
image_name = 'validation/tyoki/20180801_211605459033.jpg'

# モデルのロード
model = models.load_model(data_dir + h5_name)
model.summary()

# 画像データの準備
img = PIL.Image.open(data_dir + image_name)
img = img.resize((150, 150))
x = image.img_to_array(img)
x = x.astype('float32')/255.0
x = np.expand_dims(x, axis=0)

# 答え合わせ
hand_name = ['グー', '手が写っていない', 'その他', 'パー', 'チョキ' ]
result = model.predict(x)
your_hand_index = np.argmax(result[0])
print('あなたの出した手は... 「' + hand_name[your_hand_index] + '!」')
index = 0
for r in result[0]:
    print(hand_name[index], ' : ', r)
    index = index + 1

# 活性化を返すモデルを新たに作成
layer_outputs = [layer.output for layer in model.layers[:8]]
activation_model = models.Model(inputs=model.input, outputs=layer_outputs)
activations = activation_model.predict(x)

# 各層の名前を取得
layer_names = []
for layer in model.layers[:8]:
    layer_names.append(layer.name)

images_per_row = 8

# 特徴マップを表示
for layer_name, layer_activation in zip(layer_names, activations):
    # 特徴マップに含まれている特徴量の数
    n_features = layer_activation.shape[-1]
    # 特徴マップの形状
    size = layer_activation.shape[1]

    # この行列で活性化のチャネルをタイル表示
    n_cols = n_features // images_per_row
    display_grid = np.zeros((size * n_cols, images_per_row * size))

    # 各フィルタを1つの大きな水平グリッドでタイル表示
    for col in range(n_cols):
        for row in range(images_per_row):
            channel_image = layer_activation[0, :, :, col * images_per_row + row]
            channel_image -= channel_image.mean()
            channel_image /= channel_image.std()
            channel_image *= 64
            channel_image += 128
            channel_image = np.clip(channel_image, 0, 255).astype('uint8')
            display_grid[col * size : (col + 1) * size, row * size : (row + 1) * size] = channel_image

    # グリッドを表示
    scale = 1. / size
    plt.figure(figsize=(scale * display_grid.shape[1], scale * display_grid.shape[0]))
    plt.title(layer_name)
    plt.grid(False)
    plt.imshow(display_grid, aspect='auto', cmap='viridis')

plt.show()

ext: py

reason: import numpy


text:
import hashlib
import numpy as np
import astropy.units as u

from astropy.io import fits
from future.utils import with_metaclass
from scipy.interpolate import RegularGridInterpolator as GridInterpolate
from astropy.coordinates import SkyCoord, ICRS, BaseCoordinateFrame
from astromodels.functions.function import Function3D, FunctionMeta
from astromodels.utils.angular_distance import angular_distance_fast
from astromodels.utils.logging import setup_logger

log = setup_logger(__name__)

class GalPropTemplate_3D(with_metaclass(FunctionMeta, Function3D)):
    r"""
        description :
            Use a 3D template that has morphology and flux information.
            GalProp, DRAGON or a similar model in fits format would work.
            Only parameter is a normalization factor.
        latex : $ K $
        parameters :
            K :
                desc : normalization
                initial value : 1.0
                fix : yes
            hash :

                desc : hash of model map [needed for memoization]
                initial value : 1
                fix : yes

            lon0 :

                desc: Longitude of the center of the source
                initial value : 0.0
                min : 0.0
                max : 360.0

            lat0 :

                desc : Latitude of the center of the source
                initial value : 0.0
                min : -90.0
                max : 90.0

    """

    #__metaclass__ = FunctionMeta

    def _set_units(self, x_unit, y_unit, z_unit, w_unit):

        self.lon0.unit = x_unit
        self.lat0.unit = y_unit

        #self.K.unit = (u.MeV * u.cm**2 * u.s * u.sr) ** (-1)
        self.K.unit  = (u.sr)**(-1) #if templates are normalized, use this as the unit

    def _setup(self):

        self._frame = 'ICRS' #ICRS()
        self._map = None
        self._fitsfile = None
        self._interpmap = None

    def set_frame(self, new_frame):
        """
        Set a new frame for the coordinates (the default is ICRS J2000)
        :param new_frame: a coordinate frame from astropy
        :return: (none)
        """
        assert isinstance(new_frame, BaseCoordinateFrame)

        self._frame = new_frame

    def load_file(self, fitsfile, phi1, phi2, theta1, theta2, galactic=False, log_interp=False, ihdu=0):

        if fitsfile is None:

            log.error( "Need to specify a fits file with a template map." )

            raise RuntimeError()

        self._fitsfile = fitsfile

        p1,p2,t1,t2 = self.define_region(phi1,phi2,theta1,theta2,galactic)

        self.ramin  = p1
        self.ramax  = p2
        self.decmin = t1
        self.decmax = t2

        with fits.open(self._fitsfile) as f:

            self._delLon    = f[ihdu].header['CDELT1']
            self._delLat    = f[ihdu].header['CDELT2']
            self._delEn     = 0.2 #0.2 #f[ihdu].header['CDELT3']
            self._refLon    = f[ihdu].header['CRVAL1']
            self._refLat    = f[ihdu].header['CRVAL2']
            self._refEn     = 5 #np.log10(f[ihdu].header['CRVAL3']) # values in log10
            self._map       = f[ihdu].data
            self._refLonPix = f[ihdu].header['CRPIX1'] # reference pixel Longitude
            self._refLatPix = f[ihdu].header['CRPIX2'] # reference pixel Latitude

            self._nl        = f[ihdu].header['NAXIS1'] #longitude
            self._nb        = f[ihdu].header['NAXIS2'] #latitude
            self._ne        = f[ihdu].header['NAXIS3'] #energy

            #Create the function for the interpolation
            self._L = np.linspace(self._refLon - (self._refLonPix) * self._delLon,
                                  self._refLon + (self._nl-self._refLonPix-1) * self._delLon,
                                  self._nl)

            self._B = np.linspace(self._refLat - (self._refLatPix) * self._delLat,
                                  self._refLat + (self._nb-self._refLatPix-1) * self._delLat,
                                  self._nb)

            self._E = np.linspace(self._refEn, self._refEn+(self._ne-1)*self._delEn, self._ne)

            for i in range(len(self._E)):
                self._map[i] = (np.fliplr(self._map[i]))

            #print(np.log10(self._E))
            #print("self._E looks like: {0}".format(self._E))
            #for i in xrange(len(self._E)):
            #    self._map[i] = self._map[i].to(self.K.unit)/(np.power(10,self._E[i])*np.power(10,self._E[i])) # Map units in Mev / cm^2 s sr, changing to 1 / MeV cm^2 s sr
                #self._map[i] = (np.fliplr(self._map[i]))

            #interpolation is carried using the log scale by default as it produces better results
            # this can be disabled by changing the parameter log_interp to False
            if log_interp:

                self._F = GridInterpolate((self._E,self._B,self._L), np.log10(self._map),
                                            method="linear", bounds_error=False, fill_value=0.0)

                self._is_log10 = True

            else:

                self._F = GridInterpolate((self._E, self._B,self._L), self._map,
                                            method="linear", bounds_error=False,fill_value=0.0)

                self._is_log10 = False

            h = hashlib.sha224()
            h.update( self._map )
            self.hash = int(h.hexdigest(), 16)

    def to_dict(self, minimal=False):

        data = super(Function3D, self).to_dict(minimal)

        if not minimal:

            data['extra_setup'] = {
                                    "_fitsfile":self._fitsfile,
                                    "_frame": self._frame,
                                    "ramin": self.ramin,
                                    "ramax": self.ramax,
                                    "decmin": self.decmin,
                                    "decmax":self.decmax
                                  }

        return data

    def which_model_file(self):

        return self._fitsfile

    def evaluate(self, x, y, z, K, hash, lon0, lat0):

        if self._map is None:

            self.load_file(self._fitsfile, self.ramin, self.ramax,
                        self.decmin, self.decmax, False, True, ihdu=0)

        #transform energy from keV to MeV
        # galprop likes MeV, 3ML likes keV
        convert_val = np.log10((u.MeV.to('keV')/u.keV).value)

        # Interpolated values can be cached since we are fitting the constant K
        if self._interpmap is None:

            # We assume x and y are R.A. and Dec
            #_coord = SkyCoord(ra=x, dec=y, frame=self._frame, unit="deg")
            #b = _coord.transform_to('galactic').b.value
            #l = _coord.transform_to('galactic').l.value

            lon = x
            lat = y

            #self.lon0.value = self._refLon
            #self.lat0.value = self._refLat
            angsep = angular_distance_fast(lon0, lat0, lon, lat)

            # if only one energy is passed, make sure we can iterate just once
            if not isinstance(z, np.ndarray):

                z=np.array(z)

            energy = np.log10(z) - convert_val

            if lon.size != lat.size:

                log.error("Lon and Lat should be the same size")

                raise AttributeError()

            f_new=np.zeros([energy.size, lat.size]) # Nicola: The image is a one dimensional array

            print("Building Interpolation map....")

            num = lat.size

            for i,e in enumerate(energy):
                #print("Doing Energy E={0}".format(e))

                engs = np.repeat(e, num)

                slice_points = tuple((engs, lat, lon))

                if self._is_log10:

                    #if interpolation of flux is carried using log10 scale
                    #ensure values outside function's range remain zero, changing to linear scale
                    #makes this values evaluate to 1.
                    #i.e. if function evaluates to zero -> 10**(0) = 1. This makes the fit fail in 3ML.

                    log_slice = self._F(slice_points)

                    flux_slice = np.array(
                        [0. if x==0. else np.power(10.,x) for x in log_slice]
                    )

                else:

                    flux_slice = self._F(slice_points)

                f_new[i] = flux_slice

            log.info("Finished building interpolation map!")

            assert np.all(np.isfinite(f_new)), ("some interpolated values are wrong")

            self._interpmap=f_new

        #A = np.multiply(K,self._interpmap/(10**convert_val)) #(1000 is to change from MeV to KeV)
        A = np.multiply(K, self._interpmap) #(1000 is to change from MeV to KeV)

        return A.T

    def define_region(self, a, b, c, d, galactic=False):

        if galactic:

            lmin = a
            lmax = b
            bmin = c
            bmax = d

            _coord = SkyCoord(
                l=[lmin,lmin,lmax,lmax],
                b=[bmin, bmax, bmax, bmin],
                frame='galactic',
                unit='deg'
            )

            ramin = min(_coord.transform_to('icrs').ra.value)
            ramax = max(_coord.transform_to('icrs').ra.value)
            decmin = min(_coord.transform_to('icrs').dec.value)
            decmax = max(_coord.transform_to('icrs').dec.value)

        else:

            ramin = a
            ramax = b
            decmin = c
            decmax = d

        return ramin,ramax,decmin,decmax

    def get_boundaries(self):

        min_longitude = self.ramin
        max_longitude = self.ramax
        min_latitude = self.decmin
        max_latitude = self.decmax

        return (min_longitude, max_longitude), (min_latitude, max_latitude)

reason: import numpy,from scipy,import astropy,from astropy


text:
#!/usr/bin/env python

import numpy as np
from scipy.io import loadmat
from numpy import savez_compressed

def _mean(a, axis=None, *args, **kwargs):
    """Equivalent to numpy.mean except that the axis along which the mean is taken is not removed."""

    rv = np.mean(a, axis=axis, *args, **kwargs)

    if axis is not None:
        rv = np.expand_dims(rv, axis)

    return rv

def centre_indices(ndim=2,apron=8):
    """Returns the centre indices for the correct number of dimension
    """
    return tuple([slice(apron,-apron) for i in range(ndim)])

def summarise_mat(M, apron=8):
    """HACK to provide a 'summary' matrix consisting of the corners of the
    matrix and summed versions of the sub matrices.

    N.B. Keep this in sync with matlab/verif_m_to_npz.py.

    """
    centre = M[apron:-apron,apron:-apron,...]
    centre_sum = _mean(_mean(centre, axis=0), axis=1)

    return np.vstack((
        np.hstack((M[:apron,:apron,...], _mean(M[:apron,apron:-apron,...], axis=1), M[:apron,-apron:,...])),
        np.hstack((_mean(M[apron:-apron,:apron,...], axis=0), centre_sum, _mean(M[apron:-apron,-apron:,...], axis=0))),
        np.hstack((M[-apron:,:apron,...], _mean(M[-apron:,apron:-apron,...], axis=1), M[-apron:,-apron:,...])),
    ))

def summarise_cube(M, apron=4):
    """Provide a summary cube, extending  summarise_mat to 3D
    """
    return np.dstack(
        [summarise_mat(M[:,:,i,...], apron) for i in range(M.shape[-2])]
    )

verif_temp = loadmat('verification.mat')
verif = dict((k,v) for k, v in verif_temp.items() if (not k.startswith('_') and not k.startswith('qbgn')))
verif_cube = dict((k,v) for k, v in verif_temp.items() if (not k.startswith('_') and k.startswith('qbgn')))
del verif_temp

for idx, v in enumerate(verif['mandrill_Yh']):
    verif['mandrill_Yh_{0}'.format(idx)] = v[0]
del verif['mandrill_Yh']

for idx, v in enumerate(verif['mandrill_Yscale']):
    verif['mandrill_Yscale_{0}'.format(idx)] = v[0]
del verif['mandrill_Yscale']

for idx, v in enumerate(verif['mandrill_Yhb']):
    verif['mandrill_Yhb_{0}'.format(idx)] = v[0]
del verif['mandrill_Yhb']

for idx, v in enumerate(verif['mandrill_Yscaleb']):
    verif['mandrill_Yscaleb_{0}'.format(idx)] = v[0]
del verif['mandrill_Yscaleb']

for idx, v in enumerate(verif_cube['qbgn_Yh']):
    verif_cube['qbgn_Yh_{0}'.format(idx)] = v[0]
del verif_cube['qbgn_Yh']

for idx, v in enumerate(verif_cube['qbgn_Yscale']):
    verif_cube['qbgn_Yscale_{0}'.format(idx)] = v[0]
del verif_cube['qbgn_Yscale']

summaries = dict((k, summarise_mat(v)) for k, v in verif.items())
for k,v in verif_cube.items():
    summaries[k] = summarise_cube(v)

savez_compressed('../tests/verification.npz', **summaries)

# Convert qbgn.mat -> qbgn.npz
savez_compressed('../tests/qbgn.npz', **loadmat('qbgn.mat'))

reason: import numpy,from numpy,from scipy


text:
from numpy import *
from scipy.integrate import trapz, cumtrapz
import matplotlib
import matplotlib.pyplot as plt
import os
import sys
sys.path.append("/Users/rdodds/bin")

#
# ----------------------------------------------------------------------------
#
#    plot_set_up -- call before any plt.plot(...)
#
# ----------------------------------------------------------------------------
#
def plot_set_up( type ):
   global plot_colors
#
#      type
#       0   1 plot landscape
#       1   1 plot portrait
#
   paper_width = 11.0   # inches
   paper_height = 8.5
   params = {'legend.fontsize': 12,
      'legend.linewidth': 0.75,
      'legend.frameon': True,
      'legend.numpoints': 1,
      'figure.figsize': (paper_width,paper_height),
      'axes.linewidth': 1.125,
      'axes.titlesize': 20,       # plot title
      'axes.labelsize': 16,
      'axes.labelcolor': 'k',
      'xtick.major.size': 10,     # major tick size in points
      'xtick.minor.size': 5,      # minor tick size in points
      'xtick.major.pad': 6,       # distance to major tick label in points
      'xtick.minor.pad': 4,       # distance to the minor tick label in points
      'xtick.color': 'k',         # color of the tick labels
      'xtick.labelsize': 13,      # fontsize of the tick labels
      'ytick.major.size': 10,     # major tick size in points
      'ytick.minor.size': 5,      # minor tick size in points
      'ytick.major.pad': 6,       # distance to major tick label in points
      'ytick.minor.pad': 4,       # distance to the minor tick label in points
      'ytick.color': 'k',         # color of the tick labels
      'ytick.labelsize': 13 }     # fontsize of the tick labels
   plt.rcParams.update(params)
   plt.subplots_adjust(left=0.2, right=0.8,
                    bottom=0.2, top=0.8) # nice margins on page
   plot_colors = ["r", "g", "b", "m", "c", "k" ]
   return
#
# ----------------------------------------------------------------------------
#
#    plot_start
#
# ----------------------------------------------------------------------------
#
def plot_start():
   global plot_colors

#           1 plot on 11.0 x 8.5 (landscape)
#           has nice margins around, plot is
#           centered on page with rectangular shape
#           saved as pdf with transparent background
#
#           PPT on OS X. drag-drop .pdf onto slide.
#           scaling then changes uniformly the size of
#           all fonts, line thicknesses etc.
#
   plot_set_up( 0 )
   plt.figure()
   plot_set_up( 0 )
#
   return

#
# ----------------------------------------------------------------------------
#
#    plot_finish
#
# ----------------------------------------------------------------------------
#
def plot_finish():
   global plot_colors
#
   plt.annotate(' ',xy=(0.02, 0.92), xycoords='axes fraction')
   plt.grid(False)
#
   if os.path.isfile( plot_file_name ):
      os.remove( plot_file_name )
   plt.savefig( plot_file_name, transparent=True)
   return
#

def read_KJ_values( fname, column ):
#
 global KJ_values
 f = open(fname, "r")
 count = 0

 while True:

  while True:
   line = f.readline()
   if not line:
     break
   if line.find( "domain   KI plane stress   ") > 0:
       line = f.readline()
       line = f.readline()
       words = line.split()
       count += 1
       KJ_values[count-1,column] = float( words[2] )
       KJ_values[count-1,column] = KJ_values[count-1,column]/sqrt(1000.0)
       print KJ_values[count-1,column]

  if not line:
    break

 f.close()
 print "Normal EOF on:", filename
 return count

# ----------------------------------------------------------------------------
#
#    begin main
#
# ----------------------------------------------------------------------------
#
set_printoptions(threshold='nan')

KJ_values = zeros([1000,10], float64 )
temps = zeros(1000,float64)

filename = "./output_temp_dependent"
num_values = read_KJ_values( filename, 0 )
print "number of values: ", num_values
filename = "./output_room_temp"
num_values = read_KJ_values( filename, 1 )
print "number of values: ", num_values
filename = "./output_steady_temp"
num_values = read_KJ_values( filename, 2 )
print "number of values: ", num_values
filename = "./output_CTE_clad_zero"
num_values = read_KJ_values( filename, 3 )
print "number of values: ", num_values

now_temp = 288.0
delta_temp = -4.0
temps = zeros(1000,float64)

for j in range(num_values):
      temps[j] = now_temp
      now_temp += delta_temp
      print temps[j]

plot_start()
#
plt.ylim( [0.0, 120.0] )
plt.xlim( [288.0, 90.0] )
plt.xlabel(r'$\mathrm{Inner\ Surface\ Temperature\ (C)}$')
plt.ylabel(r'$K_J\ \ \mathrm{(MPa\sqrt{m})}$')
plt.title( " " )

plt.plot( temps[0:num_values], KJ_values[0:num_values,0], color=plot_colors[0], label="Temperature Dependent" )
plt.plot( temps[0:num_values], KJ_values[0:num_values,1], color=plot_colors[1], label="20C" )
plt.plot( temps[0:num_values], KJ_values[0:num_values,2], color=plot_colors[2], label="250C" )
plt.plot( temps[0:num_values], KJ_values[0:num_values,3], color=plot_colors[3], label="Clad CTE=0" )
#plt.plot( times_matlab, matlab, 's',markersize=3, color=plot_colors[4],label="MatLab" )

plt.legend(loc='upper left')
#
plot_file_name = 'plot.pdf'
plot_finish()

exit(0)

reason: from numpy,from scipy


text:
from collections import Counter
from typing import List, Tuple

import tensorflow as tf
import numpy as np

from tensorface.const import EMBEDDING_SIZE, UNKNOWN_CLASS

class KNN:
    def __init__(self, K=5, dist_threshold=14):
        '''
        Why such dist_threshold value?
        See notebook: notebooks/experiments_with_classification.ipynb
        :param K:
        :param dist_threshold:
        '''

        # current training data
        self.X_train = None
        self.y_train = None
        self.idx_to_lbl = None
        self.lbl_to_idx = None
        self.y_train_idx = None

        # main params
        self.dist_threshold_value = dist_threshold
        self.K = K

        # placeholders
        self.xtr = tf.placeholder(tf.float32, [None, EMBEDDING_SIZE], name='X_train')
        self.ytr = tf.placeholder(tf.float32, [None], name='y_train')
        self.xte = tf.placeholder(tf.float32, [EMBEDDING_SIZE], name='x_test')
        self.dist_threshold = tf.placeholder(tf.float32, shape=(), name="dist_threshold")

        ############ build model ############

        # model
        distance = tf.reduce_sum(tf.abs(tf.subtract(self.xtr, self.xte)), axis=1)
        values, indices = tf.nn.top_k(tf.negative(distance), k=self.K, sorted=False)
        nn_dist = tf.negative(values)
        self.valid_nn_num = tf.reduce_sum(tf.cast(nn_dist < self.dist_threshold, tf.float32))
        nn = []
        for i in range(self.K):
            nn.append(self.ytr[indices[i]])  # taking the result indexes

        # saving list in tensor variable
        nearest_neighbors = nn
        # this will return the unique neighbors the count will return the most common's index
        self.y, idx, self.count = tf.unique_with_counts(nearest_neighbors)
        self.pred = tf.slice(self.y, begin=[tf.argmax(self.count, 0)], size=tf.constant([1], dtype=tf.int64))[0]

    def predict(self, X) -> List[Tuple[str, float, List[str], List[float]]]:
        if self.X_train is None:
            # theres nothing we can do than just mark all faces as unknown...
            return [(UNKNOWN_CLASS, None, None, None) for _ in range(X.shape[0])]

        result = []
        if self.X_train is not None and self.X_train.shape[0] > 0:
            with tf.Session() as sess:
                for i in range(X.shape[0]):
                    _valid_nn_num, _pred, _lbl_idx, _counts = sess.run(
                        [self.valid_nn_num, self.pred, self.y, self.count],
                        feed_dict={
                            self.xtr: self.X_train,
                            self.ytr: self.y_train_idx,
                            self.xte: X[i, :],
                            self.dist_threshold: self.dist_threshold_value})

                    if _valid_nn_num == self.K:
                        s = _counts.sum()
                        c_lbl = []
                        c_prob = []
                        prob = None
                        for i, c in zip(_lbl_idx, _counts):
                            c_lbl.append(self.idx_to_lbl[i])
                            c_prob.append(float(c/s))
                            if _pred == i:
                                prob = float(c/s)

                        result.append((
                            self.idx_to_lbl[int(_pred)],
                            float(prob),
                            c_lbl,
                            c_prob
                        ))
                    else:
                        result.append((UNKNOWN_CLASS, None, None, None))

        return result

    def update_training(self, train_X, train_y):
        self.X_train = np.array(train_X)
        self.y_train = train_y
        self.idx_to_lbl = dict(enumerate(set(train_y)))
        self.lbl_to_idx = {v: k for k, v in self.idx_to_lbl.items()}
        self.y_train_idx = [self.lbl_to_idx[l] for l in self.y_train]

def init():
    global X, y, model
    X = []
    y = []
    model = KNN()

init()

def add(new_X, new_y):
    global X, y, model
    X.extend(new_X)
    y.extend(new_y)
    model.update_training(X, y)

def predict(X):
    global model
    return model.predict(X)

def training_data_info():
    global y
    return Counter(y)

max_stars_repo_name: CliffAgius/tf-face-recognition

max_issues_repo_name: CliffAgius/tf-face-recognition

max_forks_repo_name: CliffAgius/tf-face-recognition


text:
"""A module with platform-specific extended precision `numpy.number` subclasses.

The subclasses are defined here (instead of ``__init__.pyi``) such
that they can be imported conditionally via the numpy's mypy plugin.
"""

from typing import TYPE_CHECKING

import numpy as np
from . import (
    _80Bit,
    _96Bit,
    _128Bit,
    _256Bit,
)

if TYPE_CHECKING:
    uint128 = np.unsignedinteger[_128Bit]
    uint256 = np.unsignedinteger[_256Bit]
    int128 = np.signedinteger[_128Bit]
    int256 = np.signedinteger[_256Bit]
    float80 = np.floating[_80Bit]
    float96 = np.floating[_96Bit]
    float128 = np.floating[_128Bit]
    float256 = np.floating[_256Bit]
    complex160 = np.complexfloating[_80Bit, _80Bit]
    complex192 = np.complexfloating[_96Bit, _96Bit]
    complex256 = np.complexfloating[_128Bit, _128Bit]
    complex512 = np.complexfloating[_256Bit, _256Bit]
else:
    uint128 = NotImplemented
    uint256 = NotImplemented
    int128 = NotImplemented
    int256 = NotImplemented
    float80 = NotImplemented
    float96 = NotImplemented
    float128 = NotImplemented
    float256 = NotImplemented
    complex160 = NotImplemented
    complex192 = NotImplemented
    complex256 = NotImplemented
    complex512 = NotImplemented

max_stars_repo_name: deepyaman/numpy


text:
import os
import shutil
import tempfile

import numpy as np
import pytest

from PIL import Image

from keras_preprocessing.image import image_data_generator

@pytest.fixture(scope='module')
def all_test_images():
    img_w = img_h = 20
    rgb_images = []
    rgba_images = []
    gray_images = []
    gray_images_16bit = []
    gray_images_32bit = []
    for n in range(8):
        bias = np.random.rand(img_w, img_h, 1) * 64
        variance = np.random.rand(img_w, img_h, 1) * (255 - 64)
        # RGB
        imarray = np.random.rand(img_w, img_h, 3) * variance + bias
        im = Image.fromarray(imarray.astype('uint8')).convert('RGB')
        rgb_images.append(im)
        # RGBA
        imarray = np.random.rand(img_w, img_h, 4) * variance + bias
        im = Image.fromarray(imarray.astype('uint8')).convert('RGBA')
        rgba_images.append(im)
        # 8-bit grayscale
        imarray = np.random.rand(img_w, img_h, 1) * variance + bias
        im = Image.fromarray(imarray.astype('uint8').squeeze()).convert('L')
        gray_images.append(im)
        # 16-bit grayscale
        imarray = np.array(
            np.random.randint(-2147483648, 2147483647, (img_w, img_h))
        )
        im = Image.fromarray(imarray.astype('uint16'))
        gray_images_16bit.append(im)
        # 32-bit grayscale
        im = Image.fromarray(imarray.astype('uint32'))
        gray_images_32bit.append(im)

    return [rgb_images, rgba_images,
            gray_images, gray_images_16bit, gray_images_32bit]

def test_directory_iterator(all_test_images, tmpdir):
    num_classes = 2

    # create folders and subfolders
    paths = []
    for cl in range(num_classes):
        class_directory = 'class-{}'.format(cl)
        classpaths = [
            class_directory,
            os.path.join(class_directory, 'subfolder-1'),
            os.path.join(class_directory, 'subfolder-2'),
            os.path.join(class_directory, 'subfolder-1', 'sub-subfolder')
        ]
        for path in classpaths:
            tmpdir.join(path).mkdir()
        paths.append(classpaths)

    # save the images in the paths
    count = 0
    filenames = []
    for test_images in all_test_images:
        for im in test_images:
            # rotate image class
            im_class = count % num_classes
            # rotate subfolders
            classpaths = paths[im_class]
            filename = os.path.join(
                classpaths[count % len(classpaths)],
                'image-{}.png'.format(count))
            filenames.append(filename)
            im.save(str(tmpdir / filename))
            count += 1

    # create iterator
    generator = image_data_generator.ImageDataGenerator()
    dir_iterator = generator.flow_from_directory(str(tmpdir))

    # check number of classes and images
    assert len(dir_iterator.class_indices) == num_classes
    assert len(dir_iterator.classes) == count
    assert set(dir_iterator.filenames) == set(filenames)

    # Test invalid use cases
    with pytest.raises(ValueError):
        generator.flow_from_directory(str(tmpdir), color_mode='cmyk')
    with pytest.raises(ValueError):
        generator.flow_from_directory(str(tmpdir), class_mode='output')

    def preprocessing_function(x):
        """This will fail if not provided by a Numpy array.
        Note: This is made to enforce backward compatibility.
        """

        assert x.shape == (26, 26, 3)
        assert type(x) is np.ndarray

        return np.zeros_like(x)

    # Test usage as Sequence
    generator = image_data_generator.ImageDataGenerator(
        preprocessing_function=preprocessing_function)
    dir_seq = generator.flow_from_directory(str(tmpdir),
                                            target_size=(26, 26),
                                            color_mode='rgb',
                                            batch_size=3,
                                            class_mode='categorical')
    assert len(dir_seq) == np.ceil(count / 3.)
    x1, y1 = dir_seq[1]
    assert x1.shape == (3, 26, 26, 3)
    assert y1.shape == (3, num_classes)
    x1, y1 = dir_seq[5]
    assert (x1 == 0).all()

    with pytest.raises(ValueError):
        x1, y1 = dir_seq[14]  # there are 40 images and batch size is 3

def test_directory_iterator_class_mode_input(all_test_images, tmpdir):
    tmpdir.join('class-1').mkdir()

    # save the images in the paths
    count = 0
    for test_images in all_test_images:
        for im in test_images:
            filename = str(
                tmpdir / 'class-1' / 'image-{}.png'.format(count))
            im.save(filename)
            count += 1

    # create iterator
    generator = image_data_generator.ImageDataGenerator()
    dir_iterator = generator.flow_from_directory(str(tmpdir),
                                                 class_mode='input')
    batch = next(dir_iterator)

    # check if input and output have the same shape
    assert(batch[0].shape == batch[1].shape)
    # check if the input and output images are not the same numpy array
    input_img = batch[0][0]
    output_img = batch[1][0]
    output_img[0][0][0] += 1
    assert(input_img[0][0][0] != output_img[0][0][0])

@pytest.mark.parametrize('validation_split,num_training', [
    (0.25, 30),
    (0.50, 20),
    (0.75, 10),
])
def test_directory_iterator_with_validation_split(all_test_images,
                                                  validation_split,
                                                  num_training):
    num_classes = 2
    tmp_folder = tempfile.mkdtemp(prefix='test_images')

    # create folders and subfolders
    paths = []
    for cl in range(num_classes):
        class_directory = 'class-{}'.format(cl)
        classpaths = [
            class_directory,
            os.path.join(class_directory, 'subfolder-1'),
            os.path.join(class_directory, 'subfolder-2'),
            os.path.join(class_directory, 'subfolder-1', 'sub-subfolder')
        ]
        for path in classpaths:
            os.mkdir(os.path.join(tmp_folder, path))
        paths.append(classpaths)

    # save the images in the paths
    count = 0
    filenames = []
    for test_images in all_test_images:
        for im in test_images:
            # rotate image class
            im_class = count % num_classes
            # rotate subfolders
            classpaths = paths[im_class]
            filename = os.path.join(
                classpaths[count % len(classpaths)],
                'image-{}.png'.format(count))
            filenames.append(filename)
            im.save(os.path.join(tmp_folder, filename))
            count += 1

    # create iterator
    generator = image_data_generator.ImageDataGenerator(
        validation_split=validation_split
    )

    with pytest.raises(ValueError):
        generator.flow_from_directory(tmp_folder, subset='foo')

    train_iterator = generator.flow_from_directory(tmp_folder,
                                                   subset='training')
    assert train_iterator.samples == num_training

    valid_iterator = generator.flow_from_directory(tmp_folder,
                                                   subset='validation')
    assert valid_iterator.samples == count - num_training

    # check number of classes and images
    assert len(train_iterator.class_indices) == num_classes
    assert len(train_iterator.classes) == num_training
    assert len(set(train_iterator.filenames) &
               set(filenames)) == num_training

    shutil.rmtree(tmp_folder)

if __name__ == '__main__':
    pytest.main([__file__])

max_stars_repo_name: naitik-dodia/keras-preprocessing

max_issues_repo_name: naitik-dodia/keras-preprocessing

max_forks_repo_name: naitik-dodia/keras-preprocessing


text:
import tensorflow as tf
from slim.nets import vgg
import tensorflow.contrib.slim as slim
import numpy as np
import scipy as sp
import scipy.misc

import configargparse
import os
import sys
import datetime

# Calculates the feature correlation matrix for a 3D tensor
# (no batch dimension)
def feat_gram_matx_3d(layer, name='feat_gram_matx_3d'):
    return tf.tensordot(layer, layer, axes=[[0,1], [0,1]], name=name)

# Computes layer style losses given two feature maps from two layers
# Note that both layers should be 3D tensors (i.e. no batch dimensions)
def layer_style_loss(original_layer, generated_layer, name='layer_style_loss'):
    a = feat_gram_matx_3d(original_layer, name=name+'-original-tdot')
    g = feat_gram_matx_3d(generated_layer, name=name+'-generated-tdot')

    base_shape = tf.shape(original_layer) # should equal tf.shape(generated_layer) too
    m = tf.to_float(base_shape[0] * base_shape[1])
    n = tf.to_float(base_shape[2])

    e_unscaled = tf.reduce_sum(tf.square(g - a))

    return tf.reciprocal(4.0 * tf.square(n) * tf.square(m)) * e_unscaled

# Computes total style losses given a list of layer activation Tensors for each image
# (still no batch dimensions for the tensors)
def style_loss(original_layers, generated_layers, loss_weighting, name='style_loss'):
    layer_losses = []
    for a_layer, x_layer, w, i in zip(original_layers, generated_layers, loss_weighting, range(len(original_layers))):
        err = layer_style_loss(a_layer, x_layer, name=name+'-L{}'.format(i))
        layer_losses.append(w * err)

    return tf.reduce_sum(layer_losses, name=name+'-sum-weighted')

# Extracts the content loss by computing squared-error loss between activations
# on one layer
def content_loss(original_layer, generated_layer):
    return tf.reduce_sum(tf.square(generated_layer - original_layer) / 2.0)

# Get layer Tensors from the VGG-19 output collection dictionary
# Additionally, squeeze them all (to remove the batch dimension)
def get_style_layers(output_dict, scope_name='vgg_19'):
    #print(str(output_dict))
    return [
        tf.squeeze(output_dict[scope_name+'/conv1/conv1_1']),
        tf.squeeze(output_dict[scope_name+'/conv2/conv2_1']),
        tf.squeeze(output_dict[scope_name+'/conv3/conv3_1']),
        tf.squeeze(output_dict[scope_name+'/conv4/conv4_1']),
        tf.squeeze(output_dict[scope_name+'/conv5/conv5_1']),
    ]

# Ditto for the content representation layer
def get_content_layer(output_dict, scope_name='vgg_19'):
    return tf.squeeze(output_dict[scope_name+'/conv4/conv4_2'])

def model_preprocessing(input_img_data):
    ret = tf.to_float(input_img_data)
    return ret - tf.constant([123.68, 116.779, 103.939]) # Subtract mean pixel values

def decode_image_pipeline(input_filedata, input_filename):
    fnparts = os.path.splitext(input_filename)
    ext = fnparts[1]
    if ext == '.png':
        return tf.image.decode_png(input_filedata, channels=3)
    elif ext == '.jpg' or ext == '.jpeg':
        return tf.image.decode_jpeg(input_filedata, channels=3)

# Performs precomputation for content and style images.
def build_graph_stage1(args):
    print('Building precompute graph...')
    sys.stdout.flush()

    file_contents = {}

    filenames = [args.content_image]
    file_weights = [] # don't include content weighting
    for fp in args.style_image:
        pair = fp.split(':')
        filenames.append(pair[0])
        file_weights.append(float(pair[1]))

    image_tensors = []
    for filename in filenames:
        handle = open(filename, 'rb')
        data = handle.read()
        handle.close()

        fd_tensor = tf.constant(data, dtype=tf.string, name='data-'+os.path.basename(filename))
        decoded_tensor = decode_image_pipeline(fd_tensor, filename)
        model_input = model_preprocessing(decoded_tensor)
        model_input = tf.reshape(model_input, [args.image_size[0], args.image_size[1], 3])

        image_tensors.append(model_input)

    stacked_images = tf.parallel_stack(image_tensors)

    _, activations = vgg.vgg_19(stacked_images, is_training=False, spatial_squeeze=False)

    content_layer = get_content_layer(activations)[0]   # 1st batch element == content image
    style_layers = get_style_layers(activations)        # All other batch elements == style images

    n_style_images = len(filenames)-1
    style_activations = [[] for i in range(n_style_images)] # 2D list of activations for each style image

    for batched_layer in style_layers:
        unbatched = tf.unstack(batched_layer, num=len(filenames))[1:]
        for j, activation in enumerate(unbatched):
            style_activations[j].append(activation)

    model_vars = tf.get_collection(tf.GraphKeys.GLOBAL_VARIABLES, scope='vgg_19')
    saver = tf.train.Saver(model_vars)

    print("Launching precompute graph...")
    sys.stdout.flush()

    with tf.Session() as stage1:
        print('Loading model weights from {} and initializing global variables...'.format(args.vgg_19_checkpoint))
        sys.stdout.flush()

        stage1.run(tf.global_variables_initializer())
        saver.restore(stage1, args.vgg_19_checkpoint)

        print('Precomputing content and style image activations...')
        sys.stdout.flush()

        fetches = stage1.run([content_layer, style_activations])

        actual_content_layer = fetches[0]
        actual_style_activations = fetches[1]

    return actual_content_layer, actual_style_activations, file_weights

# Builds the portions of the graph that actually carry out the main computation
def build_graph_stage2(args, content_layer, style_activations, file_weights):
    generated_image_data = tf.placeholder(tf.float32, shape=args.image_size[0]*args.image_size[1]*3, name='generated_image')
    generated_image = tf.reshape(generated_image_data, [1]+args.image_size+[3])

    img_summary = tf.summary.image('Generated Image', generated_image, collections=[])

    generated_image = model_preprocessing(generated_image)

    # Run the generated image through the ConvNet model...
    _, generated_outputs = vgg.vgg_19(generated_image, is_training=False, spatial_squeeze=False)

    # Extract important layers from the model collection dictionary.
    # Since we build stage2 in its own graph, we use the default scope name for finding the layers.
    generated_content_layer = get_content_layer(generated_outputs, 'vgg_19')
    generated_style_layers = get_style_layers(generated_outputs, 'vgg_19')

    # Convert the precomputed layers to constant tensors
    constant_content_layer = tf.squeeze(tf.constant(content_layer, dtype=tf.float32, name='constant_content_input'))
    constant_style_layers = []
    for i, style_image in enumerate(style_activations):
        t = []
        for j, layer in enumerate(style_image):
            t.append(
                tf.squeeze(tf.constant(layer, dtype=tf.float32, name='constant_style_activation-{}-{}'.format(i,j)))
            )
        constant_style_layers.append(t)

    # Calculate loss subcomponents w/ weights
    wcl = args.content_loss_weight * content_loss(constant_content_layer, generated_content_layer)
    tf.summary.scalar('Weighted Content Loss', wcl)

    weighted_losses = [wcl]
    for i, image_layers in enumerate(constant_style_layers):
        image_style_loss = style_loss(image_layers, generated_style_layers, [0.2] * 5, name='style-loss-I{}'.format(i))
        weighted_loss = image_style_loss * file_weights[i]

        tf.summary.scalar('Weighted Style Loss {}'.format(i), weighted_loss)
        weighted_losses.append(weighted_loss)

    # Compute overall loss
    total_loss = tf.reduce_sum(weighted_losses)
    tf.summary.scalar('Total Loss', total_loss)

    # Compute gradients
    gen_grad = tf.gradients(total_loss, generated_image)
    flattened_gradients = tf.reshape(gen_grad, [-1])
    tf.summary.scalar('Mean Gradient', tf.reduce_mean(gen_grad))

    # Auxillary stuff
    model_vars = tf.get_collection(tf.GraphKeys.GLOBAL_VARIABLES, scope='vgg_19')
    saver = tf.train.Saver(model_vars)

    increment_step = tf.assign_add(tf.contrib.framework.get_or_create_global_step(), 1)

    outputs = {
        'loss': total_loss,
        'grad': flattened_gradients,
        'summary': tf.summary.merge_all(),
        'image_summary': img_summary,
        'increment_step': increment_step,
    }

    return generated_image_data, outputs, saver

def init_compute(args, cluster=None, server=None):
    precompute_graph = tf.Graph()

    dev_fn = None
    session_target = None
    if args.distributed:
        dev_fn = tf.train.replica_device_setter(
            worker_device="/job:worker/task:{:d}".format(args.task_index),
            cluster=cluster
        )
        session_target = server.target

    with slim.arg_scope([slim.variable, slim.model_variable], device='/cpu:0'):
        with precompute_graph.as_default():
            content_layer, style_acts, style_weights = build_graph_stage1(args)

    print('Building main compute graph.')
    sys.stdout.flush()

    with slim.arg_scope([slim.variable, slim.model_variable], device='/cpu:0'):
        with tf.device(dev_fn):
            image_input_tensor, outputs, model_saver = build_graph_stage2(args, content_layer, style_acts, style_weights)

    #tf.train.write_graph(precompute_graph, args.logdir, 'precompute-graph.pbtxt')
    #tf.train.write_graph(tf.get_default_graph(), args.logdir, 'compute-graph.pbtxt')

    print("Launching main graph...")
    sys.stdout.flush()

    sess = tf.Session(target=session_target)

    summ_writer = tf.summary.FileWriter(args.logdir, graph=sess.graph)

    print("Initializing variables...")
    sys.stdout.flush()
    sess.run(tf.global_variables_initializer())

    print('Reloading model weights from {}...'.format(args.vgg_19_checkpoint))
    sys.stdout.flush()
    model_saver.restore(sess, args.vgg_19_checkpoint)

    # Prepare initial image and bounds for L-BFGS-B
    img = None
    if args.continue_optimization is not None:
        img = scipy.misc.imread(args.continue_optimization).flatten()
        print("Loaded interim image! Continuing...")
        sys.stdout.flush()
    else:
        img = np.random.uniform(0.0, 255.0, size=args.image_size[0]*args.image_size[1]*3)
    bounds = np.zeros((args.image_size[0]*args.image_size[1]*3, 2))
    bounds[:] = (0.0, 255.0)

    # Get our loss function for L-BFGS-B:
    current_iteration = 1

    ###
    def get_loss(generated_img):
        nonlocal current_iteration

        if current_iteration % args.iteration_output_freq == 0:
            print("[{}] - beginning iteration {}".format(str(datetime.datetime.now()), current_iteration))
            sys.stdout.flush()

        loss, grads, summ, step = sess.run(
            (outputs['loss'], outputs['grad'], outputs['summary'], outputs['increment_step']),
            feed_dict={ image_input_tensor: generated_img }
        )

        current_iteration = step+1
        summ_writer.add_summary(summ, step)

        if args.interim_file_template is not None and args.interim_save_freq is not None:
            if step % args.interim_save_freq == 0:
                interm_img = np.reshape(generated_img, args.image_size + [3])
                scipy.misc.imsave(args.interim_file_template + '{}.png'.format(step), interm_img)

                img_summ = outputs['image_summary'].eval(session=sess, feed_dict={ image_input_tensor: generated_img })
                summ_writer.add_summary(img_summ, step)

        if step % args.iteration_output_freq == 0:
            print('[{}] Iteration {} - loss: {}'.format(str(datetime.datetime.now()), step, loss))
            sys.stdout.flush()

        return np.float64(loss), np.float64(grads) # SciPy doesn't like it when you pass float32's to minimize...
    ###

    return sess, img, bounds, get_loss

def main_compute_loop(args):
    ps_hosts = None
    cluster = None
    server = None

    if not os.path.exists(args.logdir):
        if os.path.exists(os.path.dirname(args.logdir)):
            os.makedirs(args.logdir)
        else:
            raise FileNotFoundError("Logging base directory {} does not exist.".format(os.path.dirname(args.logdir)))

    if args.interim_file_template is not None:
        interim_img_dir = os.path.dirname(args.interim_file_template)
        if not os.path.exists(interim_img_dir):
            if os.path.exists(os.path.dirname(interim_img_dir)):
                os.makedirs(interim_img_dir)
            else:
                raise FileNotFoundError("Intermediate image base directory {} does not exist.".format(os.path.dirname(interim_img_dir)))

    if args.distributed:
        print("Initializing distributed mode.")
        sys.stdout.flush()

        ps_hosts = args.ps_hosts.split(',')
        cluster = tf.train.ClusterSpec({
            "ps": ps_hosts,
            "worker": [args.worker_host],
        })
        server = tf.train.Server(cluster, job_name=args.job, task_index=args.task_index)

    if args.distributed and args.job == 'ps':
        print("Starting parameter server.")
        sys.stdout.flush()
        server.join()
    else:
        sess, init_img, bounds, loss_fn = init_compute(args, cluster, server)

        print("Beginning image generation.")
        sys.stdout.flush()

        final_img, final_loss, info_dict = sp.optimize.fmin_l_bfgs_b(
            loss_fn,
            init_img,
            bounds=bounds,
            iprint=args.iteration_output_freq,
            maxiter=args.max_iterations,
        )

        print("Saving output...")
        sys.stdout.flush()

        output_img = np.reshape(final_img, args.image_size + [3])
        scipy.misc.imsave(args.outfile, output_img)

# Main script code below:
parser = configargparse.ArgumentParser(description='Performs neural style transfer.', default_config_files=['~/.nst-settings'])
parser.add_argument('-c', '--config', is_config_file=True, help='Configuration file path')

parser.add_argument('--image_size', required=True, type=int, nargs=2, help='Height and width of images used for content, style, and output (i.e. --image_size 512 512)')
parser.add_argument('--content_image', required=True, help='Input content image')
parser.add_argument('--outfile', required=True, help='File to write completed image to (image dimensions will be determined from content image)')

parser.add_argument('--style_image', action='append', help='An input style image path:weight pair (i.e. my_image.jpg:1.0). This option can be specified multiple times.')
parser.add_argument('--content_loss_weight', default=8e-4, type=float, help='Alpha parameter for loss calculation')

parser.add_argument('--distributed', action='store_true', help='Run in distributed mode')
parser.add_argument('--ps_hosts', help='List of parameter server hostname:port pairs, comma-separated.')
parser.add_argument('--worker_host', help='The worker node\'s hostname:port pair')
parser.add_argument('--job', help='Job name for this instance.')
parser.add_argument('--task_index', type=int, default=0, help='Task index for specified job.')

parser.add_argument('--continue_optimization', default=None, help='Path to an intermediate image to resume from (begin from random data if not set)')
parser.add_argument('--vgg_19_checkpoint', default='vgg_19.ckpt', help='Checkpoint file containing VGG-19 model weights')
parser.add_argument('--logdir', default='nst-logdir', help='Directory to save summaries to')
parser.add_argument('--interim_file_template', default='nst-intermediate-', help='Intermediate file template (iteration number and .png extension will be appended)')
parser.add_argument('--interim_save_freq', default=10, type=int, help='How often to save intermediate images during optimization')
parser.add_argument('--iteration_output_freq', default=1, type=int, help='How often to print information to console during optimization')
parser.add_argument('--max_iterations', default=15000, type=int, help='Maximum iterations to run for output')

args = parser.parse_args()
main_compute_loop(args)

max_stars_repo_name: stmobo/Machine-Learning

max_issues_repo_name: stmobo/Machine-Learning

max_forks_repo_name: stmobo/Machine-Learning

reason: import numpy,import scipy


text:
import numpy as np
import math
import matplotlib.pyplot as plt
from PIL import Image, ImageDraw
import torch
from torch.nn.functional import softmax
from torchvision import transforms
import timm

def visualize_filters(model, filter_name = None, max_filters = 64, size = 128, figsize = (16, 16), save_path = None):
    name, weights = next(model.named_parameters())

    for layer_name, layer_weights in model.named_parameters():
        if layer_name == filter_name:
            name = layer_name
            weights = layer_weights

    w_size = weights.size()
    merged_weights = weights.reshape(w_size[0] * w_size[1], w_size[2], w_size[2]).detach().numpy()
    out_chs = merged_weights.shape[0]

    if out_chs > max_filters:
        merged_weights = merged_weights[torch.randperm(out_chs)[:max_filters]]
        out_chs = max_filters

    sqrt = int(math.sqrt(out_chs))
    fig, axs = plt.subplots(sqrt, sqrt, figsize = figsize)

    if not size:
        size = merged_weights.shape[2]

    for i in range(sqrt ** 2):
        weight = merged_weights[i]
        scaled = scale(weight)
        resized = transforms.Resize((size, size))(Image.fromarray(scaled))
        plot_idx = int(i / sqrt), i % sqrt

        axs[plot_idx].imshow(resized, cmap = 'gray')
        axs[plot_idx].set_yticks([])
        axs[plot_idx].set_xticks([])

    if save_path:
        fig.savefig(save_path)

def visualize_activations(model, module, img_path, max_acts = 64, figsize = (16, 16), save_path = None):
    img_t = preprocess_image(img_path)
    acts = [0]

    def hook_fn(self, input, output):
        acts[0] = output

    handle = module.register_forward_hook(hook_fn)
    out = model(img_t)
    handle.remove()
    acts = acts[0][0].detach().numpy()

    if acts.shape[0] > max_acts:
        acts = acts[torch.randperm(acts.shape[0])[:max_acts]]

    sqrt = int(math.sqrt(acts.shape[0]))
    fig, axs = plt.subplots(sqrt, sqrt, figsize = figsize)

    for i in range(sqrt ** 2):
        scaled = scale(acts[i])

        plt_idx = int(i / sqrt), i % sqrt
        axs[plt_idx].imshow(scaled, cmap = 'gray')
        axs[plt_idx].set_yticks([])
        axs[plt_idx].set_xticks([])

    if save_path:
        fig.savefig(save_path)

def maximally_activated_patches(model, img_path, patch_size = 448, stride = 100, num_patches = 5, figsize = (16, 16),
                                    device = 'cuda', save_path = None):
    model.eval()
    model.to(device)

    img = Image.open(img_path).convert('RGB')
    img_mean = int(np.mean(np.asarray(img)))
    img_t = preprocess_image(img_path).to(device)

    with torch.no_grad():
        out = model(img_t)

    probs = softmax(out[0], dim = 0)
    max_index = probs.argmax()
    orig = probs[max_index]

    dim1 = int(((img.size[0] - patch_size) / stride) + 1)
    dim2 = int(((img.size[1] - patch_size) / stride) + 1)
    diff = []

    for i in range(dim1 * dim2):
        img_copy = img.copy()

        x0, y0, x1, y1 = gen_coords(i, patch_size, stride, dim1, dim2)

        draw = ImageDraw.Draw(img_copy)
        draw.rectangle([x0, y0, x1, y1], fill = (img_mean, img_mean, img_mean))

        occ_img_t = preprocess_image(img_copy).to(device)

        with torch.no_grad():
            out = model(occ_img_t)

        probs_occ = softmax(out[0], dim = 0)
        diff.append(abs(orig - probs_occ[max_index].item()))

    diff = np.array(diff)
    top_indices = diff.argsort()[-num_patches:]
    fig, axs = plt.subplots(int(num_patches / 5), min(num_patches, 5), figsize = figsize)

    for i, idx in enumerate(top_indices):
        img_copy = img.copy()
        x0, y0, x1, y1 = gen_coords(idx, patch_size, stride, dim1, dim2)

        if num_patches > 5:
            plot_idx = int(i / 5), i % 5
        else:
            plot_idx = i

        axs[plot_idx].imshow(np.asarray(img_copy.crop((x0, y0, x1, y1))))
        axs[plot_idx].set_yticks([])
        axs[plot_idx].set_xticks([])

def saliency_map(model, img_path, figsize = (16, 16), device = 'cuda', save_path = None):
    model.eval()
    model.to(device)

    img = Image.open(img_path).convert('RGB')
    img_t = preprocess_image(img).to(device)
    img_np = scale(img_t.detach().cpu()[0].permute(1, 2, 0).numpy())
    img_t.requires_grad = True

    out = model(img_t)
    max_out = out.max()
    max_out.backward()
    saliency, _ = torch.max(img_t.grad.data.abs(), dim = 1)
    saliency = saliency.squeeze(0)
    saliency_img = saliency.detach().cpu().numpy()

    fig, axs = plt.subplots(1, 2, figsize = figsize)
    axs[0].imshow(img_np)
    axs[1].imshow(saliency_img, cmap = 'gray')
    axs[0].set_yticks([])
    axs[0].set_xticks([])
    axs[1].set_yticks([])
    axs[1].set_xticks([])

    if save_path:
        plt.savefig(save_path)

def generate_image(model, target_class, epochs, min_prob, lr, weight_decay, step_size = 100, gamma = 0.6,
                        noise_size = 224, p_freq = 50, init = torch.randn, device = 'cuda', figsize = (6, 6)):

    noise = init([1, 3, noise_size, noise_size]).to(device)
    noise.requires_grad = True
    model = model.to(device)
    opt = torch.optim.SGD([noise], lr = lr, weight_decay = weight_decay)
    scheduler = torch.optim.lr_scheduler.StepLR(opt, step_size, gamma = gamma)

    for i in range(1, epochs + 1):
        opt.zero_grad()
        outs = model(noise)
        p = softmax(outs[0], dim = 0)[target_class]

        if i % p_freq == 0 or i == epochs:
            print('Epoch: {} Confidence score for class {}: {}'.format(i, target_class, p))

        if p > min_prob:
            print('Reached {} confidence score in epoch {}. Stopping early.'.format(p, i))
            break

        obj = - outs[0][target_class]
        obj.backward()
        opt.step()
        scheduler.step()

    fig, axs = plt.subplots(1, figsize = figsize)
    img_np = postprocess_image(noise)
    axs.imshow(img_np)
    axs.set_xticks([])
    axs.set_yticks([])

    return noise

def fool_model(model, img_path, target_class, epochs, min_prob, lr,
                        step_size, gamma, p_freq = 50, init = torch.randn, device = 'cuda', figsize = (6, 6)):

    orig_img = preprocess_image(img_path).to(device)
    orig_img.requires_grad = True
    model = model.to(device)
    opt = torch.optim.SGD([orig_img], lr = lr)
    scheduler = torch.optim.lr_scheduler.StepLR(opt, step_size, gamma = gamma)

    for i in range(1, epochs + 1):
        opt.zero_grad()
        outs = model(orig_img)
        p = softmax(outs[0], dim = 0)[target_class]

        if i % p_freq == 0 or i == epochs:
            print('Epoch: {} Confidence score for class {}: {}'.format(i, target_class, p))

        if p > min_prob:
            print('Reached {} confidence score in epoch {}. Stopping early.'.format(p, i))
            break

        obj = - outs[0][target_class]
        obj.backward()
        opt.step()
        scheduler.step()

    fig, axs = plt.subplots(1, figsize = figsize)
    img_np = postprocess_image(orig_img)
    axs.imshow(img_np)
    axs.set_xticks([])
    axs.set_yticks([])

    return orig_img

def feature_inversion_helper(module, orig_img, epochs, lr, step_size, gamma, mu, noise_size = 224, init = torch.randn,
                                device = 'cuda'):
    acts = [0]
    def hook_fn(self, input, output):
        acts[0] = output

    handle = module.register_forward_hook(hook_fn)
    _ = model(orig_img)
    orig_features = acts[0]

    noise = init([1, 3, noise_size, noise_size]).to(device)
    noise.requires_grad = True
    opt = torch.optim.SGD([noise], lr = lr)
    scheduler = torch.optim.lr_scheduler.StepLR(opt, step_size, gamma = gamma)

    for i in range(epochs):
        opt.zero_grad()
        _ = model(noise)
        curr = acts[0]

        loss = ((orig_features - curr) ** 2).sum() + mu * total_variation_regularizer(noise)
        loss.backward(retain_graph = True)

        opt.step()
        scheduler.step()

    handle.remove()
    return noise

def feature_inversion(model, module_dict, img_path, epochs, lr, step_size = 100, gamma = 0.6, mu = 1e-1,
                          device = 'cuda', figsize = (16, 16)):

    orig_img = preprocess_image(img_path).to(device)
    model = model.to(device)
    model.eval()
    recreated_imgs = []

    for module in modules:
        recreated_imgs.append(feature_inversion_helper(module, orig_img, epochs = epochs,
                                lr = lr, step_size = step_size, gamma = gamma, mu = mu, device = device))

    fig, axs = plt.subplots(1, len(recreated_imgs), figsize = figsize)

    for i in range(len(recreated_imgs)):
        axs[i].imshow(postprocess_image(recreated_imgs[i]))

def deep_dream(model, module, img_path, epochs, lr, step_size = 100, gamma = 0.6, device = 'cuda', figsize = (12, 12)):

    img_t = preprocess_image(img_path).to(device)
    img_t.requires_grad = True
    opt = torch.optim.SGD([img_t], lr = lr)
    scheduler = torch.optim.lr_scheduler.StepLR(opt, step_size, gamma = gamma)
    acts = [0]

    def hook_fn(self, input, output):
        acts[0] = output

    model.to(device)
    model.eval()
    handle = module.register_forward_hook(hook_fn)

    for i in range(epochs):
        opt.zero_grad()
        _ = model(img_t)
        loss = -acts[0].norm()
        loss.backward()
        opt.step()
        scheduler.step()

    handle.remove()

    fig, axs = plt.subplots(1, figsize = figsize)
    img_np = postprocess_image(img_t)
    axs.imshow(img_np)
    axs.set_xticks([])
    axs.set_yticks([])

    return img_t


text:
import math
import os
import time
import webbrowser
import wmi
from ctypes import cast, POINTER
import cv2
import numpy as np
from comtypes import CLSCTX_ALL
from cvzone.HandTrackingModule import HandDetector
from pycaw.pycaw import AudioUtilities, IAudioEndpointVolume

# volume metrics
devices = AudioUtilities.GetSpeakers()
interface = devices.Activate(IAudioEndpointVolume._iid_, CLSCTX_ALL, None)
volume = cast(interface, POINTER(IAudioEndpointVolume))
volRange = volume.GetVolumeRange()
minVol = volRange[0]
maxVol = volRange[1]
vol = 0
volBar = 400
volPer = 0

# brightness metrics
minBrightness = 0
maxBrightness = 100

class HandDetectorClass():
    def __init__(self):
        super(HandDetectorClass, self).__init__()
        self.cap = cv2.VideoCapture(0)
        self.detector = HandDetector(detectionCon=0.8, maxHands=2)
        self.x = 0

    def readVideoFromCamera(self):
        while True:
            self.success, self.img = self.cap.read()
            self.hands, self.img = self.detector.findHands(self.img)
            if self.hands:
                self.InitializeHand()
                try:
                    self.DoActivity(self.x)
                except Exception as e:
                    print(e)
            cv2.imshow("Image", self.img)
            cv2.waitKey(1)

    def InitializeHand(self):
        self.hand1 = self.hands[0]
        self.lmList1 = self.hand1["lmList"]
        self.bbox1 = self.hand1["bbox"]
        self.centerPoint1 = self.hand1["center"]
        self.handType1 = self.hand1["type"]
        self.fingers1 = self.detector.fingersUp(self.hand1)

        if len(self.hands) == 2:
            self.hand2 = self.hands[1]
            self.lmList2 = self.hand2["lmList"]
            self.bbox2 = self.hand2["bbox"]
            self.centerPoint2 = self.hand2["center"]
            self.handType2 = self.hand2["type"]
            self.fingers2 = self.detector.fingersUp(self.hand2)
            try:
                self.length, self.info, self.img = self.detector.findDistance(self.centerPoint1, self.centerPoint2,
                                                                              self.img)
                self.x = int(self.length)
            except Exception as e:
                print("-----------")

        elif len(self.hands) == 1 and self.handType1 == "Left":
            self.VolumeControls(self.lmList1)

        elif len(self.hands) == 1 and self.handType1 == "Right":
            self.BrightnessControls(self.lmList1)

    def VolumeControls(self, lmList1):
        x1, y1 = self.lmList1[4][1], self.lmList1[4][2]
        x2, y2 = self.lmList1[8][1], self.lmList1[8][2]
        cx, cy = (x1 + x2) // 2, (y1 + y2) // 2

        cv2.circle(self.img, (x1, y1), 15, (255, 0, 255), cv2.FILLED)
        cv2.circle(self.img, (x2, y2), 15, (255, 0, 255), cv2.FILLED)
        cv2.line(self.img, (x1, y1), (x2, y2), (255, 0, 255), 3)
        cv2.circle(self.img, (cx, cy), 15, (255, 0, 255), cv2.FILLED)

        length = math.hypot(x2 - x1, y2 - y1)
        # print(length)
        # Hand range 50 - 300
        # Volume Range -65 - 0
        vol = np.interp(length, [50, 300], [minVol, maxVol])
        volBar = np.interp(length, [50, 300], [400, 150])
        volPer = np.interp(length, [50, 300], [0, 100])
        # print(int(length), vol)
        volume.SetMasterVolumeLevel(vol, None)

        if length < 50:
            cv2.circle(self.img, (cx, cy), 15, (0, 255, 0), cv2.FILLED)

        cv2.rectangle(self.img, (50, 150), (85, 400), (255, 0, 0), 3)
        cv2.rectangle(self.img, (50, int(volBar)), (85, 400), (255, 0, 0), cv2.FILLED)
        cv2.putText(self.img, f'{int(volPer)} %', (40, 450), cv2.FONT_ITALIC, 1, (255, 0, 0), 3)

    def BrightnessControls(self, lmList1):
        x11, y11 = self.lmList1[4][1], self.lmList1[4][2]
        x22, y22 = self.lmList1[8][1], self.lmList1[8][2]
        cx1, cy1 = (x11 + x22) // 2, (y11 + y22) // 2
        cv2.circle(self.img, (x11, y11), 15, (255, 0, 255), cv2.FILLED)
        cv2.circle(self.img, (x22, y22), 15, (255, 0, 255), cv2.FILLED)
        cv2.line(self.img, (x11, y11), (x22, y22), (255, 0, 255), 3)
        cv2.circle(self.img, (cx1, cy1), 15, (255, 0, 255), cv2.FILLED)
        length = math.hypot(x22 - x11, y22 - y11)
        brightness = np.interp(length, [50, 300], [minBrightness, maxBrightness])
        c = wmi.WMI(namespace='wmi')
        methods = c.WmiMonitorBrightnessMethods()[0]
        methods.WmiSetBrightness(brightness, 0)
        # print(int(brightness))
        BrightnessBar = np.interp(length, [50, 300], [400, 150])
        cv2.rectangle(self.img, (50, 150), (85, 400), (255, 0, 0), 3)
        cv2.rectangle(self.img, (50, int(BrightnessBar)), (85, 400), (255, 0, 0), cv2.FILLED)
        cv2.putText(self.img, f'{int(brightness)} %', (40, 450), cv2.FONT_ITALIC, 1, (255, 0, 0), 3)

    def DoActivity(self, x):
        if self.x in range(400, 500, 2):
            self.openWhatsApp()

    def openWhatsApp(self):
        webbrowser.open("https://web.whatsapp.com")

obj = HandDetectorClass()
obj.readVideoFromCamera()

max_stars_repo_name: RakiRam/PC-Automation-Using-Gestures

max_issues_repo_name: RakiRam/PC-Automation-Using-Gestures

max_forks_repo_name: RakiRam/PC-Automation-Using-Gestures


text:
import numpy as np
import pandas as pd
from math import acos, atan2
import shutil
import pathlib
import os
import re
from glob import glob
from datetime import datetime
from random import random
import warnings
from Bio.PDB import PDBParser, HSExposureCB
from Bio.PDB.DSSP import DSSP
from itertools import groupby, chain, combinations
import operator
from copy import deepcopy
from collections import ChainMap

__location__ = os.path.realpath(os.path.join(os.getcwd(), os.path.dirname(__file__)))
forster_radius = 54.0  # forster radius in Angstrom for Cy3-Cy5, according to Murphy et al. 2004

# --- IO ---
def parse_output_dir(out_dir, clean=False):
    out_dir = os.path.abspath(out_dir) + '/'
    if clean:
        shutil.rmtree(out_dir, ignore_errors=True)
    pathlib.Path(out_dir).mkdir(parents=True, exist_ok=True)
    return out_dir

def parse_input_dir(in_dir, pattern=None, regex=None):
    if type(in_dir) != list: in_dir = [in_dir]
    out_list = []
    for ind in in_dir:
        if not os.path.exists(ind):
            raise ValueError(f'{ind} does not exist')
        if os.path.isdir(ind):
            ind = os.path.abspath(ind)
            if pattern is not None: out_list.extend(glob(f'{ind}/**/{pattern}', recursive=True))
            else: out_list.extend(glob(f'{ind}/**/*', recursive=True))
        else:
            if pattern is None: out_list.append(ind)
            elif pattern.strip('*') in ind: out_list.append(ind)
    if regex is not None:
        out_list = [fn for fn in out_list if re.search(regex, fn)]
    return out_list

def print_timestamp():
    return datetime.now().strftime('%y-%m-%d_%H:%M:%S')

def get_angle(v1, v2):
    # pv = np.cross(v1,v2)
    # pv = pv / np.linalg.norm(pv)
    #
    # pv2 = np.cross(v1, pv)
    # pv2 = pv2 / np.linalg.norm(pv2)
    #
    # # new basis
    # nb = np.vstack((v1, pv, pv2))
    #
    # v2_nb = np.dot(nb, v2)
    #
    # return np.degrees(atan2(np.linalg.det([v1,v2,pv]), np.dot(v1, v2)))
    return np.degrees(acos(max(-1, min(1, np.dot(v1, v2) / (np.linalg.norm(v1) * np.linalg.norm(v2))))))

def get_abs_dihedral(ca1, ca2, cb1, cb2):
    """
    Get absolute (!) dihedral from coordinates
    source: https://stackoverflow.com/questions/20305272/dihedral-torsion-angle-from-four-points-in-cartesian-coordinates-in-python
    """
    p0 = cb1
    p1 = ca1
    p2 = ca2
    p3 = cb2

    b0 = -1.0*(p1 - p0)
    b1 = p2 - p1
    b2 = p3 - p2

    # normalize b1 so that it does not influence magnitude of vector
    # rejections that come next
    b1 = b1 / np.linalg.norm(b1)

    # vector rejections
    # v = projection of b0 onto plane perpendicular to b1
    #   = b0 minus component that aligns with b1
    # w = projection of b2 onto plane perpendicular to b1
    #   = b2 minus component that aligns with b1
    v = b0 - np.dot(b0, b1)*b1
    w = b2 - np.dot(b2, b1)*b1

    # angle between v and w in a plane is the torsion angle
    # v and w may not be normalized but that's fine since tan is y/x
    x = np.dot(v, w)
    y = np.dot(np.cross(b1, v), w)
    return np.abs(np.degrees(np.arctan2(y, x)))

def get_compound_properties(pdb_id, pdb_fn, target_residue_list):
    """
    Retrieve selected properties from first chain (alphabetically) in a pdb file, in order:
    - compound name
    - synonyms
    - length of chain
    - list of number of residues present of given type(s)
    """
    pdb_parser = PDBParser()
    with warnings.catch_warnings():
        warnings.simplefilter('ignore')
        pdb_obj = pdb_parser.get_structure(pdb_id, pdb_fn)
    chain = sorted(list(pdb_obj.get_chains()), key=lambda x: x.get_id())[0]
    chain_id = chain.get_id().lower()
    header = [ch for ch in pdb_obj.header['compound'].values() if chain_id in ch['chain']][0]
    aa_str = ''.join([aa_dict_31.get(res.get_resname(), '') for res in chain.get_residues()])
    res_occurances = [aa_str.count(res) for res in target_residue_list]
    return header['molecule'], header['synonym'], len(chain), res_occurances

def get_tagged_resi(mm, tagged_resn, aa_sequence):
    tagged_resi = {}
    tagged_resi_list = [0]  # N-terminus is always labeled
    for target_resn in tagged_resn:
        resi_list = []
        cur_mm = mm.get(target_resn, {target_resn: 1.0})
        for resi, resn in enumerate(aa_sequence):
            if resi not in tagged_resi_list and random() < cur_mm.get(resn, 0.0):
                resi_list.append(resi)
        tagged_resi[target_resn] = resi_list
        tagged_resi_list.extend(resi_list)
    return tagged_resi

# --- additional higher-dim numpy functions ---
def np2d_argmin(mat):
    d1 = mat.argmin()
    return d1 // mat.shape[1], d1 % mat.shape[1]

def inNd(a, b, assume_unique=False):
    if a.ndim == 1:
        a = np.expand_dims(a, 0)
    if b.ndim == 1:
        b = np.expand_dims(b, 0)
    if a.shape[1] != b.shape[1]:
        # return None
        raise ValueError(f'dim 1 of a ({a.shape[1]}) and b ({b.shape[1]}) do not match')
    # a = np.asarray(a, order='C')
    # b = np.asarray(b, order='C')
    a = a.ravel().view((np.void, a.dtype.itemsize*a.shape[1]))
    b = b.ravel().view((np.void, b.dtype.itemsize*b.shape[1]))
    return np.in1d(a, b, assume_unique)

def contains_double_coords(coords):
    sorted_coords = coords[np.lexsort(coords.T),:]
    double_bool = np.invert(np.any(np.diff(sorted_coords, axis=0), 1))
    return np.any(double_bool)

# --- lattice structure mutation functions ---
def pick_random_coords(coords):
    if coords.size == 0:
        raise ValueError('coordinates array is empty')
    elif coords.shape[0] == 1:
        return coords[0]
    else:
        return coords[np.random.randint(coords.shape[0]), :]

def get_neighbors(c, d=1):
    neighbors = np.tile(c, (6, 1))
    neighbors += np.row_stack((np.eye(3, dtype=int) * d, np.eye(3, dtype=int) * -1) * d)
    return neighbors

def get_diagonal_neighbors(c, d=1):
    neighbors = np.tile(c, (12, 1))
    m1 = np.abs(np.eye(3, dtype=int) - 1) * d
    m2 = np.array([[0, 1, -1],[1,0,-1],[1,-1,0]], dtype=int) * d
    m3 = m1 * -1
    m4 = m2 * -1
    mods = np.vstack((m1, m2, m3, m4))
    neighbors[:, :-1] += mods
    return neighbors

# def get_rotmat(dim):
#     th = 0.5 * np.pi * (1, -1)[dim < 0]
#     d = abs(dim)
#     if d == 1:
#         return np.array([[1, 0, 0, 0],
#                         [0, np.cos(th), -1*np.sin(th), 0],
#                         [0, np.sin(th), np.cos(th), 0],
#                         [0, 0, 0, 1]], dtype=int)
#     elif d == 2:
#         return np.array([[np.cos(th), 0, np.sin(th), 0],
#                         [0, 1, 0, 0],
#                         [-1*np.sin(th), 0, np.cos(th), 0],
#                         [0, 0, 0, 1]], dtype=int)
#     elif d == 3:
#         return np.array([[np.cos(th), -1*np.sin(th), 0, 0],
#                         [np.sin(th), np.cos(th), 0, 0],
#                         [0, 0, 1, 0],
#                         [0, 0, 0, 1]], dtype=int)

# def get_transmat(c, dir):
#     em = np.eye(4, dtype=int)
#     em[:-1, 3] = c[:-1] * dir
#     return em

def pdb_coord(c):
    c_str_list = []
    for nc in c:
        c_str = f'{nc:.3f}'
        c_str = ' ' * (8 - len(c_str)) + c_str
        c_str_list.append(c_str)
    return ''.join(c_str_list)

def get_FRET_efficiency(dist):
    """
    return FRET efficiency for a given distance, between cy3-cy5
    """
    return 1.0/((dist / forster_radius) ** 6 + 1)

def get_FRET_distance(efret):
    return (efret ** -1 - 1) ** (1/6) * forster_radius

def get_pairs_mat(path):
    pairs_mat = pd.read_csv(path, sep='\s+', index_col=0, header=None)
    if len(pairs_mat.index[0]) == 3:
        pairs_mat.index = [aa_dict_31.get(str(aa).upper(), str(aa).upper()) for aa in pairs_mat.index]
    pairs_mat.columns = pairs_mat.index
    pairs_mat.sort_index(inplace=True)
    pairs_mat = pairs_mat.reindex(sorted(pairs_mat.columns), axis=1)
    pairs_mat.loc[:, 'TAG'] = 0.0
    pairs_mat.loc['TAG', :] = 0.0
    return pairs_mat

# --- residue name conversion dictionaries ---
aa_dict = {
    'A': 'ALA',
    'C': 'CYS',
    'D': 'ASP',
    'E': 'GLU',
    'F': 'PHE',
    'G': 'GLY',
    'H': 'HIS',
    'I': 'ILE',
    'K': 'LYS',
    'L': 'LEU',
    'M': 'MET',
    'N': 'ASN',
    'P': 'PRO',
    'Q': 'GLN',
    'R': 'ARG',
    'S': 'SER',
    'T': 'THR',
    'V': 'VAL',
    'W': 'TRP',
    'Y': 'TYR',
    'X': 'TAG'
}

# For coloring html plots of structures:
# - blue: pos
# - red: neg
# - grey: polar
# - yellow: hydrophobic
aa_to_hp_dict = {
    'A': 'yellow',
    'C': 'grey',
    'D': 'red',
    'E': 'red',
    'F': 'yellow',
    'G': 'yellow',
    'H': 'blue',
    'I': 'yellow',
    'K': 'blue',
    'L': 'yellow',
    'M': 'yellow',
    'N': 'grey',
    'P': 'yellow',
    'Q': 'grey',
    'R': 'blue',
    'S': 'grey',
    'T': 'grey',
    'V': 'yellow',
    'W': 'yellow',
    'Y': 'grey'
}

aa_dict_31 = {v: k for k, v in aa_dict.items()}

atm_mass_dict = {
    'C': 12.0107,
    'O': 15.9994,
    'N': 14.0067,
    'S': 32.065
}

aa_mass_dict = {  # todo: where did I get these?? Does not quite match other sources...
    'A': 71.00779,
    'C': 103.1429,
    'D': 115.0874,
    'E': 129.11398,
    'F': 147.17386,
    'G': 57.05132,
    'H': 137.13928,
    'I': 113.15764,
    'K': 128.17228,
    'L': 113.15764,
    'M': 131.19606,
    'N': 114.10264,
    'P': 97.11518,
    'Q': 128.12922,
    'R': 156.18568,
    'S': 87.0773,
    'T': 101.10388,
    'V': 99.13106,
    'W': 186.2099,
    'Y': 163.17326
}

# cm_dist_df = pd.read_pickle(f'{__location__}/data/aa_pairwise_dist_df.pkl')

rotmat_dict = {
    1: np.array([[1,0,0],[0,0,-1],[0,1,0]]),
    -1: np.array([[1,0,0],[0,0,1],[0,1,0]]),
    2: np.array([[0,0,1],[0,1,0],[-1,0,0]]),
    -2: np.array([[0,0,1],[0,1,0],[1,0,0]]),
    3: np.array([[0,-1,0],[1,0,0],[0,0,1]]),
    -3: np.array([[0,1,0],[1,0,0],[0,0,1]])
}

# --- pair potentials functions ---

def remove_nonstd_residues(chain):
    """
    Remove non-standard residues that were in first/last position, check for nonstandard residues in other positions
    """
    res_list = list(chain.get_residues())
    caps = [res_list[0], res_list[-1]]
    while not all( [aa_dict_31.get(cap.get_resname(), False) for cap in caps]):
        [chain.detach_child(cap.get_id()) for cap in caps if not aa_dict_31.get(cap.get_resname(), False)]
        res_list = list(chain.get_residues())
        caps = [res_list[0], res_list[-1]]
    resbool_list = [ aa_dict_31.get(res.get_resname(), False) for res in  chain.get_residues()]
    return all(resbool_list)

def get_pairs_counts(pdb_fn, count_tup, parallel_bool, out_queue, disorder_only):
    try:
        pair_count_mat, ss_mat, count_dict = count_tup
        mat_update = np.eye(2, dtype=int)
        pdb_id = os.path.splitext(os.path.basename(pdb_fn))[0]
        structure = PDBParser(QUIET=True).get_structure(pdb_id, pdb_fn)
        model = structure[0]
        nb_residues = len(list(model.get_residues()))

        resn_dict = {res.id[1]: aa_dict_31[res.get_resname()] for res in model.get_residues()}

        # --- Secondary structure counts ---
        # todo: sheet df index changed!! No longer 0-based!
        helix_df, sheet_df = parse_ss_info(pdb_fn)
        helix_list, sheet_list = [], []
        if len(helix_df):
            helix_list = np.concatenate([np.arange(tup.resi_start, tup.resi_end+1) for _, tup in helix_df.iterrows()]) + 1
        if len(sheet_df):
            sheet_list = np.concatenate([np.arange(tup.resi_start, tup.resi_end+1) for _, tup in sheet_df.iterrows()]) + 1
        ss_dict = {res.id[1]: 'L' for res in model.get_residues()}
        for ssi in helix_list:
            if ssi in  ss_dict: ss_dict[ssi] = 'H'
        for ssi in sheet_list:
            if ssi in ss_dict: ss_dict[ssi] = 'S'
        for ssi in ss_dict:
            ss_mat.loc[resn_dict[ssi], ss_dict[ssi]] += 1

        # --- neighbor counts ---
        # hse = HSExposureCB(model, radius=6.6, offset=1)
        # hse

        # --- old stuff ---
        # # hse = HSExposureCB(model, radius=8, offset=1)
        # pdb_dict = {resi: {'resn': aa_dict_31[res.get_resname()], 'nb_adjacent': sum(hse[res.full_id[2:]][:2]),
        #                    'coords': [], 'ss': ss_dict.get(res.id[1], 'L')} for resi, res in enumerate(model.get_residues())}
        # cb_coords = np.empty((nb_residues, 3))
        # cb_coords[:,:] = np.nan
        #
        # # Collect c-beta coordinates and residue counts
        # for ri, res in enumerate(model.get_residues()):
        #     resname = aa_dict_31[res.get_resname()]
        #     assert pdb_dict[ri]['resn'] == resname
        #     if resname == 'G':
        #         cb = [atm for atm in res.get_atoms() if atm.get_name() == 'CA']
        #         coords = cb[0].get_coord()
        #         pdb_dict[ri]['coords'] = coords
        #     else:
        #         cb = [atm for atm in res.get_atoms() if atm.get_name() == 'CB']
        #         if not len(cb):
        #             if not out_queue: return pair_count_mat, count_dict
        #             out_queue.put((pair_count_mat, count_dict))
        #             return
        #         coords = cb[0].get_coord()
        #         cb_coords[ri, :] = coords
        #         pdb_dict[ri]['coords'] = coords
        #
        # # Update counts
        # for ri in pdb_dict:
        #     ss_mat.loc[pdb_dict[ri]['resn'], pdb_dict[ri]['ss']] += 1
        #     if disorder_only:
        #         if pdb_dict[ri]['ss'] != 'L': continue
        #     resname = pdb_dict[ri]['resn']
        #     # nb_neighbors = round(4 * (1 - pdb_dict[ri]['asa']))
        #     nb_neighbors = min(4, pdb_dict[ri]['nb_adjacent'])
        #     nb_water = 4 - nb_neighbors
        #     if nb_neighbors != 0:  # neighbor interactions
        #         exclude_idx = np.arange(max(0, ri-2), min(ri+2, nb_residues))
        #         coord_dist = np.linalg.norm(pdb_dict[ri]['coords'] - cb_coords, axis=1)
        #         proximity_order = np.argsort(coord_dist)
        #         proximity_order = proximity_order[np.invert(np.in1d(proximity_order, exclude_idx))]
        #         neighbors = [pdb_dict[ni] for ni in proximity_order[:nb_neighbors]]
        #         for neighbor in neighbors:
        #             if disorder_only:
        #                 if neighbor['ss'] != '-': continue
        #             pair_count_mat.loc[(resname, neighbor['resn']), (neighbor['resn'], resname)] += mat_update
        #     if nb_water != 0:  # water interactions
        #         count_dict['HOH'] += nb_water
        #         pair_count_mat.loc[(resname, 'HOH'), ('HOH', resname)] += mat_update * nb_water
        #     count_dict[resname] += 1
        if parallel_bool:
            out_queue.put((pair_count_mat, ss_mat, count_dict))
        else:
            return pair_count_mat, ss_mat, count_dict
    except:
        if parallel_bool:
            out_queue.put(count_tup)
        else:
            return count_tup

# --- to generate mock ss data from a pdb file, as if from ss prediction tool ---
def get_mock_ss_df(pdb_fn, resi2idx_dict, skip_outer_helix_indices=False):
    nb_res = len(resi2idx_dict)
    helix_txt_list = []
    sheet_txt_list = []
    with open(pdb_fn, 'r') as fh:
        for line in fh.readlines():
            if line.startswith('HELIX'):
                helix_txt_list.append(line)
            elif line.startswith('SHEET'):
                sheet_txt_list.append(line)
    helix_df = parse_helix_info(helix_txt_list)
    sheet_df = parse_sheet_info(sheet_txt_list)
    if skip_outer_helix_indices:
        helix_df = helix_df.query('length >=3').copy()
        helix_df.resi_start = helix_df.resi_start + 1
        helix_df.resi_end = helix_df.resi_end - 1
        helix_df.length = helix_df.length - 2

    helix_dict = {}
    for helix_id, tup in helix_df.iterrows():
        res_idx_list, helix_coords = put_helix_on_lattice(tup)
        res_idx_list = [resi2idx_dict[hidx] for hidx in res_idx_list if hidx in resi2idx_dict]
        helix_dict[np.min(res_idx_list)] = (res_idx_list, helix_coords)
    ss_df = pd.DataFrame(np.tile([0.05, 0.05, 0.90], (nb_res, 1)), columns=['H', 'S', 'L'], index=np.arange(nb_res))
    for _, tup in helix_df.iterrows():
        if tup.resi_end - tup.resi_start < 4: continue
        ah_idx = np.arange(resi2idx_dict[tup.resi_start], resi2idx_dict[tup.resi_end] - 3)
        ss_df.loc[ah_idx, :] = [0.90, 0.05, 0.05]
    for _, tup in sheet_df.iterrows():
        sh_idx = np.arange(resi2idx_dict[tup.resi_start], resi2idx_dict[tup.resi_end] + 1)
        ss_df.loc[sh_idx, :] = [0.05, 0.90, 0.05]
    ss_df = np.log(1 / ss_df - 1)
    return ss_df, helix_dict

def put_helix_on_lattice(tup):
    # generate all possible conformations for type of helix
    mods = helix_type_mod_dict[tup.type]
    nb_steps = mods[0].shape[0]
    coord_list = []
    for mod in mods:
        coords = np.zeros((tup.length, 3))
        for i in range(tup.length - 1):
            coords[i + 1, :3] = coords[i, :3] + mod[i % nb_steps, :]
        coords_mirrored = get_rotated_poses(coords)
        coord_list.extend(coords_mirrored)
    return np.arange(tup.resi_start, tup.resi_end + 1), coord_list

helix_type_mod_dict = {
    # 1: [np.array([[2, 2, 2],  # compressed
    #              [-2, 2, 0],
    #              [-2, -2, 2],
    #              [2, -2, 2]])],
    # 5: [np.array([[2, 2, 2],  # todo find better implementation for 3-10
    #              [-2, 2, 0],
    #              [-2, -2, 2],
    #              [2, -2, 2]])],
    1: [np.array([[2, -2, 2],  # uncompressed, full length
                 [2, 2, 2],
                 [-2, 2, 2],
                 [-2, -2, 2]])],
    5: [np.array([[2, -2, 2],
                 [2, 2, 2],
                 [-2, 2, 2],
                 [-2, -2, 2]])],

    # 5: [np.array([[1, 0, 0],
    #              [-1, 0, 1],
    #              [0, 1, -1]]),
    #     np.array([[1, 0, 1],
    #               [-1, 0, 0],
    #               [0, 1, -1]]),
    #     np.array([[1, 0, 1],
    #               [-1, 1, 0],
    #               [0, 0, -1]]),
    #     ]
}

def get_rotated_poses(coords):
    out_list = []
    for rot in rotmats:
        coords_rot = coords.copy()
        coords_rot[:, :3] = np.matmul(coords[:, :3], rot)
        out_list.append(coords_rot.copy())
    return out_list

rotmats_single = [
    np.array([[0, -1, 0],
              [1, 0, 0],
              [0, 0, 1]]),
    np.array([[0, 0, 1],
              [0,1,0],
              [-1,0,0]]),
    np.array([[1,0,0],
              [0,0,-1],
              [0,1,0]]),
    np.array([[1, 0, 0],
              [0, -1, 0],
              [0, 0, -1]]),
    np.array([[-1,0,0],
              [0,1,0],
              [0,0,-1]]),
    np.array([[-1,0,0],
              [0,-1,0],
              [0,0,1]]),
    np.array([[1,0,0],
              [0,0,1],
              [0,-1,0]]),
    np.array([[0,0,-1],
              [0,1,0],
              [1,0,0]]),
    np.array([[0,1,0],
              [-1,0,0],
              [0,0,1]])
]

# mirror_dims = list(combinations([0,1,2], 2)) + [tuple([i]) for i in range(3)] + [(0, 1, 2)]
# mirror_dims_dual = list(combinations([0,1,2], 2))
rotmats_2 = [np.matmul(rotmats_single[md[0]], rotmats_single[md[1]]) for md in list(combinations(np.arange(len(rotmats_single)), 2))]
rotmats_3 = [np.matmul(np.matmul(rotmats_single[md[0]], rotmats_single[md[1]]), rotmats_single[md[2]]) for md in list(combinations(np.arange(len(rotmats_single)), 3))]
rotmats_w_doubles = list(chain.from_iterable([rotmats_single, rotmats_2, rotmats_3]))
rotmats = []
for rm in rotmats_w_doubles:
    dists = [np.linalg.norm(rm - el) for el in rotmats]
    if not len(dists) or 0 not in dists:
        rotmats.append(rm)

# ---

def parse_ss_info(pdb_fn):
    helix_txt_list = []
    sheet_txt_list = []
    with open(pdb_fn, 'r') as fh:
        for line in fh.readlines():
            if line.startswith('HELIX'):
                helix_txt_list.append(line)
            elif line.startswith('SHEET'):
                sheet_txt_list.append(line)
    helix_df = parse_helix_info(helix_txt_list)
    sheet_df = parse_sheet_info(sheet_txt_list)
    return helix_df, sheet_df

def parse_sheet_info(sheet_txt_list):
    """
    Parse list of strings of PDB SHEET records to pandas data frame with 0-based residue indices
    """
    sheet_df = pd.DataFrame(
        columns=['resi_start', 'resi_end', 'orientation', 'resi_hb_cur', 'resi_hb_prev'],
    index=pd.MultiIndex(levels=[[], []], codes=[[], []], names=['strand_id', 'sheet_id']))
    for line in sheet_txt_list:
        ll = line.split()
        if ll[1] == '1':
            sheet_df.loc[(int(ll[1]),ll[2]), :] = [int(line[22:26]), int(line[33:37]), int(line[38:40]),
                                                   pd.NA, pd.NA]
        elif ll[1].isdigit():
            sheet_df.loc[(int(ll[1]),ll[2]), :] = [int(line[22:26]), int(line[33:37]), int(line[38:40]),
                                                   int(line[50:54]), int(line[65:69])]
        else:
            break
    return sheet_df

def parse_helix_info(helix_txt_list):
    """
    Parse list of strings of PDB HELIX records to a pandas data frame with 0-based residue indices
    """
    helix_df = pd.DataFrame(
        columns=['resi_start', 'resi_end', 'type', 'length']
    )
    for line in helix_txt_list:
        helix_df.loc[int(line[7:10]), :] = [int(line[21:25]), int(line[33:37]), int(line[39:40]), int(line[72:76])]
    return helix_df

def most_frequent_element(id_list):
    """
    Return most frequent list element, or with lowest index if tied
    from: https://stackoverflow.com/questions/1518522/find-the-most-common-element-in-a-list
    """
    id_list_sorted = sorted((idl, i) for i, idl in enumerate(id_list))
    groups = groupby(id_list_sorted, key=operator.itemgetter(0))
    def _auxfun(g):
        item, iterable = g
        count = 0
        min_index = len(id_list)
        for _, where in iterable:
          count += 1
          min_index = min(min_index, where)
        # print 'item %r, count %r, minind %r' % (item, count, min_index)
        return count, -min_index
    return max(groups, key=_auxfun)[0]

def generate_pdb(bb_coord_list, resnames, res_coord_list=None, first_n=None):
    # Create pdb file
    coord_txt_list = []
    conect_txt_list = []
    atm_nb = 1
    prev_ca = 1
    for ri, (bb_coord, resname) in enumerate(zip(bb_coord_list, resnames)):
        if ri == 0 and first_n is not None:
            # also register first N
            n_str = pdb_coord(first_n)
            coord_txt_list.append(
                f'HETATM{str(atm_nb).rjust(5)}  N   {resname} A{str(ri).rjust(4)}    {n_str}  1.00  1.00           N\n')
            atm_nb += 1
            conect_txt_list.append(f"CONECT{str(atm_nb - 1).rjust(5)}{str(atm_nb).rjust(5)}\n")

        ca_str = pdb_coord(bb_coord)
        coord_txt_list.append(
            f'HETATM{str(atm_nb).rjust(5)}  CA  {resname} A{str(ri).rjust(4)}    {ca_str}  1.00  1.00           C\n')
        conect_txt_list.append(
            f"CONECT{str(prev_ca).rjust(5)}{str(atm_nb).rjust(5)}\n")  # connect to previous CA (or first N)
        prev_ca = atm_nb
        atm_nb += 1
        if res_coord_list is not None and resname != 'GLY':
            cb_str = pdb_coord(res_coord_list[ri])
            coord_txt_list.append(
                f'HETATM{str(atm_nb).rjust(5)}  CB  {resname} A{str(ri).rjust(4)}    {cb_str}  1.00  1.00           C\n')
            atm_nb += 1
            conect_txt_list.append(f"CONECT{str(atm_nb - 2).rjust(5)}{str(atm_nb - 1).rjust(5)}\n")
    coord_txt_list.append(f'TER{str(atm_nb - 2).rjust(5)}      {resnames[-1]} A{str(len(resnames)).rjust(4)}\n')
    conect_txt = ''.join(conect_txt_list)
    coord_txt = ''.join(coord_txt_list)
    pdb_txt = coord_txt + conect_txt
    return pdb_txt

def get_cm(atm_list):
        return np.sum(np.row_stack([atm.get_coord() * atm.mass for atm in atm_list]), axis=0) / \
               np.sum([atm.mass for atm in atm_list])

def put_cb_on_lattice(cb_array, cb_mod_array, ca_array, ca_array_lat, resname_list, n_coord, n1_dist_unit):
    nb_res = len(resname_list)
    for ri in range(nb_res):
        if ri == 0:
            # also put first N on lattice
            neighbors_unit = get_neighbors(ca_array_lat[ri, :]).astype(int)
            neighbors = neighbors_unit * n1_dist_unit
            n_neigh_idx = np.argmin(np.linalg.norm(neighbors[:, :3] - n_coord, axis=1))
            n_coord_unit = neighbors_unit[n_neigh_idx, :]
            n_coord_lat = neighbors[n_neigh_idx, :]

        # apply same translation on CB as was done for CA
        if resname_list[ri] == 'GLY': continue
        transl = ca_array_lat[ri, :] - ca_array[ri, :]
        cb_coord = cb_array[ri, :] + transl

        # generate all possible positions for CB, pick closest to real coords
        neighbors = get_neighbors(ca_array_lat[ri, :]).astype(int)
        if ri == 0:
            adjacent_cas = np.vstack((ca_array_lat[ri+1, :], n_coord_unit))
        elif ri == nb_res - 1:
            adjacent_cas = np.expand_dims(ca_array_lat[ri-1, :], 0)
        else:
            adjacent_cas = np.vstack((ca_array_lat[ri-1, :], ca_array_lat[ri+1, :]))
        neighbors = neighbors[np.invert(inNd(neighbors[:, :3], adjacent_cas[:, :3])), :]
        cb_coord_overshot = neighbors[np.argmin(np.linalg.norm(neighbors[:, :3] - cb_coord[:3], axis=1)), :3]
        cb_mod_array[ri, :3] = cb_coord_overshot - ca_array_lat[ri, :3]

def list_sheet_series(ss_seq):
    dict_list = []
    in_sheet = False
    cur_list = []
    for ti, ss in enumerate(ss_seq):
        if ss == 'S':
            cur_list.append(ti)
        else:
            if len(cur_list): dict_list.append({cl: deepcopy(cur_list) for cl in cur_list})
            cur_list = []
    out_dict = dict(ChainMap(*dict_list))
    return out_dict