text:
import numpy as np
import pandas as pd
from scipy.stats import multivariate_normal
from statsmodels.api import GLM, families

from .bins import atleast_2d

def estimate_movement_std(position):
    '''Estimates the movement standard deviation based on position.

    WARNING: Need to use on original position, not interpolated position.

    Parameters
    ----------
    position : ndarray, shape (n_time, n_position_dim)

    Returns
    -------
    movement_std : ndarray, shape (n_position_dim,)

    '''
    position = atleast_2d(position)
    is_nan = np.any(np.isnan(position), axis=1)
    position = position[~is_nan]
    movement_std = []
    for p in position.T:
        fit = GLM(p[:-1], p[1:], family=families.Gaussian()).fit()
        movement_std.append(np.sqrt(fit.scale))
    return np.array(movement_std)

def _normalize_row_probability(x):
    '''Ensure the state transition matrix rows sum to 1
    '''
    x /= x.sum(axis=1, keepdims=True)
    x[np.isnan(x)] = 0
    return x

def _fix_zero_bins(movement_bins):
    '''If there is no data observed for a column, set everything to 1 so
    that it will have equal probability
    '''
    n_bins = movement_bins.shape[0]
    movement_bins[movement_bins.sum(axis=1) == 0] = 1 / n_bins
    return movement_bins

def empirical_movement(position, edges, is_training=None, replay_speed=20,
                       position_extent=None):
    '''Estimate the probablity of the next position based on the movement
     data, given the movment is sped up by the
     `replay_speed`

    Place cell firing during a hippocampal replay event is a "sped-up"
    version of place cell firing when the animal is actually moving.
    Here we use the animal's actual movements to constrain which place
    cell is likely to fire next.

    Parameters
    ----------
    position : ndarray, shape (n_time, n_position_dims)
    edges : sequence
        A sequence of arrays describing the bin edges along each dimension.
    is_training : None or bool ndarray, shape (n_time,), optional
    replay_speed : int, optional
        How much the movement is sped-up during a replay event
    position_extent : sequence, optional
        A sequence of `n_position_dims`, each an optional (lower, upper)
        tuple giving the outer bin edges for position.
        An entry of None in the sequence results in the minimum and maximum
        values being used for the corresponding dimension.
        The default, None, is equivalent to passing a tuple of
        `n_position_dims` None values.

    Returns
    -------
    transition_matrix : ndarray, shape (n_position_bins, n_position_bins)

    '''
    if is_training is None:
        is_training = np.ones((position.shape[0]), dtype=np.bool)
    position = atleast_2d(position)[is_training]
    movement_bins, _ = np.histogramdd(
        np.concatenate((position[1:], position[:-1]), axis=1),
        bins=edges * 2, range=position_extent)
    original_shape = movement_bins.shape
    n_position_dims = position.shape[1]
    shape_2d = np.product(original_shape[:n_position_dims])
    movement_bins = _normalize_row_probability(
        movement_bins.reshape((shape_2d, shape_2d), order='F'))
    movement_bins = np.linalg.matrix_power(movement_bins, replay_speed)

    return movement_bins

def random_walk(place_bin_centers, movement_var, is_track_interior,
                replay_speed=20):
    '''Zero mean random walk.

    This version makes the sped up gaussian without constraints and then
    imposes the constraints.

    Parameters
    ----------
    place_bin_centers : ndarray, shape (n_bins, n_position_dims)
    movement_var : float,
    is_track_interior : bool ndarray, shape (n_x_bins, n_y_bins)
    replay_speed : int

    Returns
    -------
    transition_matrix : ndarray, shape (n_bins, n_bins)

    '''
    transition_matrix = np.stack(
        [multivariate_normal(mean=bin, cov=movement_var * replay_speed).pdf(
            place_bin_centers)
         for bin in place_bin_centers], axis=1)
    is_track_interior = is_track_interior.ravel(order='F')
    transition_matrix[~is_track_interior] = 0.0
    transition_matrix[:, ~is_track_interior] = 0.0

    return _normalize_row_probability(transition_matrix)

def _center_arm(row, bin_labels, before_gaussian, after_gaussian):
    arm_ind = bin_labels.loc[(bin_labels == 'Right Arm')].index
    n_samples = min(arm_ind.size, after_gaussian.size)
    row[arm_ind[:n_samples]] = after_gaussian[:n_samples] / 2

    arm_ind = bin_labels.loc[(bin_labels == 'Left Arm')].index
    n_samples = min(arm_ind.size, after_gaussian.size)
    row[arm_ind[:n_samples]] = after_gaussian[:n_samples] / 2
    return row

def _left_arm(row, bin_labels, before_gaussian, after_gaussian):
    arm_ind = bin_labels.loc[(bin_labels == 'Center Arm')].index
    n_samples = min(arm_ind.size, before_gaussian.size)
    row[arm_ind[-n_samples:]] = before_gaussian[-n_samples:] / 2

    arm_ind = bin_labels.loc[(bin_labels == 'Right Arm')].index
    n_samples = min(arm_ind.size, before_gaussian.size)
    row[arm_ind[:n_samples]] = before_gaussian[-n_samples:][::-1] / 2
    return row

def _right_arm(row, bin_labels, before_gaussian, after_gaussian):
    arm_ind = bin_labels.loc[(bin_labels == 'Center Arm')].index
    n_samples = min(arm_ind.size, before_gaussian.size)
    row[arm_ind[-n_samples:]] = before_gaussian[-n_samples:] / 2

    arm_ind = bin_labels.loc[(bin_labels == 'Left Arm')].index
    n_samples = min(arm_ind.size, before_gaussian.size)
    row[arm_ind[:n_samples]] = before_gaussian[-n_samples:][::-1] / 2
    return row

_ARM_FUNCS = {
    'Left Arm': _left_arm,
    'Right Arm': _right_arm,
    'Center Arm': _center_arm,
}

def w_track_1D_random_walk(position, place_bin_edges, place_bin_centers,
                           labels, movement_var, is_track_interior,
                           replay_speed=200):
    position = position.squeeze()
    place_bin_edges = place_bin_edges.squeeze()
    place_bin_centers = place_bin_centers.squeeze()
    bin_inds = np.digitize(position, bins=place_bin_edges) - 1
    n_bins = place_bin_edges.size - 1
    bin_inds[bin_inds >= n_bins] = n_bins - 1
    bin_labels = pd.Series(labels).groupby(
        bin_inds).unique().apply(lambda s: s[0])

    transition_matrix = []

    for bin_ind in np.arange(n_bins):
        bin_center = place_bin_centers[bin_ind]
        row = np.zeros((n_bins,))
        try:
            bin_label = bin_labels.loc[bin_ind]
            is_same_bin = (bin_labels == bin_label)
            is_same_bin_ind = is_same_bin[is_same_bin].index
            is_same_bin = np.zeros_like(
                place_bin_centers.squeeze(), dtype=np.bool)
            is_same_bin[is_same_bin_ind] = True

            gaussian = multivariate_normal(
                mean=bin_center, cov=movement_var * replay_speed).pdf(
                place_bin_centers)
            row[is_same_bin] = gaussian[is_same_bin]
            before_gaussian = gaussian[:np.nonzero(is_same_bin)[0][0]]
            after_gaussian = gaussian[np.nonzero(is_same_bin)[0][-1]:]
            row = _ARM_FUNCS[bin_label](row, bin_labels, before_gaussian,
                                        after_gaussian)
        except KeyError:
            pass
        transition_matrix.append(row)

    transition_matrix = np.stack(transition_matrix, axis=1)
    is_track_interior = is_track_interior.ravel(order='F')
    transition_matrix[~is_track_interior] = 0.0
    transition_matrix[:, ~is_track_interior] = 0.0
    return _normalize_row_probability(transition_matrix)

def random_walk_on_track_graph(
    place_bin_centers, movement_mean, movement_var,
    place_bin_center_ind_to_node,
    distance_between_nodes
):
    state_transition = np.zeros(
        (place_bin_centers.size, place_bin_centers.size))
    gaussian = multivariate_normal(mean=movement_mean, cov=movement_var)

    for bin_ind1, node1 in enumerate(place_bin_center_ind_to_node):
        for bin_ind2, node2 in enumerate(place_bin_center_ind_to_node):
            try:
                state_transition[bin_ind1, bin_ind2] = gaussian.pdf(
                    distance_between_nodes[node1][node2]
                )
            except KeyError:
                # bins not on track interior will be -1 and not in distance
                # between nodes
                continue

    return state_transition

reason: import numpy,from scipy,from statsmodels


text:
import cv2
import numpy as np

device = 0
cap=cv2.VideoCapture(device)
k=np.ones((3,3),np.uint8)
t0=cap.read()[1]
t1=cap.read()[1]

while(True):
        d=cv2.absdiff(t1,t0)
        grey=cv2.cvtColor(d, cv2.COLOR_BGR2GRAY)
        blur = cv2.GaussianBlur(grey, (3,3), 0)
        ret, th = cv2.threshold( blur, 15, 255, cv2.THRESH_BINARY )
        dilated=cv2.dilate(th,k,iterations=2)
        # if need additional value, use and underscore at the beginning
        contours, hierarchy = cv2.findContours(
		dilated,cv2.RETR_TREE,cv2.CHAIN_APPROX_SIMPLE)
        t2=t0
        cv2.drawContours(t2, contours, -1, (0,255,0), 2 )
        cv2.imshow("Output:", t2 )
        t0=t1
        t1=cap.read()[1]

        # click on camera window and press ESC
        if cv2.waitKey(5) == 27:
                break

cap.release()
cv2.destroyAllWindows()

ext: py

max_stars_repo_name: msmcs-robotics/pi-cv-vision

max_issues_repo_name: msmcs-robotics/pi-cv-vision

max_forks_repo_name: msmcs-robotics/pi-cv-vision

reason: import numpy


text:
#!/usr/bin/env python3

"""
Implementation of some test on the vanishing point detection methods.
"""

###########
# Imports #
###########

import cv2
import numpy as np
import os
import sys

from PIL import Image

current = os.path.dirname(os.path.realpath(__file__))
parent = os.path.dirname(os.path.dirname(current))
sys.path.append(parent)

from analysis.vanishing_point import VPClassic, VPEdgelets  # noqa: E402

########
# Main #
########

def main(
    method_id: str = 'classic',
    inpt_pth: str = 'image.png',
    outpt_pth: str = 'results/'
):
    # Load method
    methods = {
        'classic': VPClassic,
        'edgelets': VPEdgelets
    }

    method = methods.get(method_id)(export=True)

    # Open image
    img = Image.open(inpt_pth)
    img = np.array(img)

    print(f'Image dimensions: {img.shape}')

    # Get results
    guess, exports = method.detect(img)

    # Export results
    os.makedirs(os.path.dirname(outpt_pth), exist_ok=True)

    for i, export in enumerate(exports):
        cv2.imwrite(
            os.path.join(outpt_pth, f'{i}.png'),
            export
        )

    # Display vanishing point guess
    print(f'Vanishing point guess: {guess}')

if __name__ == '__main__':
    import argparse

    parser = argparse.ArgumentParser(
        description='Test or evaluate vanishing point detectors.'
    )

    parser.add_argument(
        '-m',
        '--method',
        type=str,
        default='classic',
        choices=['classic', 'edgelets'],
        help='method to use'
    )

    parser.add_argument(
        '-i',
        '--input',
        type=str,
        default='image.png',
        help='path to input file'
    )

    parser.add_argument(
        '-o',
        '--output',
        type=str,
        default='results/',
        help='path to output file or folder'
    )

    args = parser.parse_args()

    main(
        method_id=args.method,
        inpt_pth=args.input,
        outpt_pth=args.output
    )

max_stars_repo_name: meurissemax/autonomous-drone

max_issues_repo_name: meurissemax/autonomous-drone

max_forks_repo_name: meurissemax/autonomous-drone


text:
# encoding: utf-8
# network file -> build basic pipline and decoder for Dynamic Network
# @author: yanwei.li
import torch
import torch.nn as nn
import torch.nn.functional as F
import numpy as np
from typing import Dict
from dl_lib.layers import Conv2d, ShapeSpec, get_norm
from dl_lib.modeling.nn_utils import weight_init
from dl_lib.modeling.postprocessing import sem_seg_postprocess
from dl_lib.structures import ImageList
from dl_lib.modeling.dynamic_arch import cal_op_flops
from dl_lib.modeling.dynamic_arch import cal_op_flops
from dl_lib.modeling.dynamic_arch.dynamic_backbone import build_dynamic_backbone
from dl_lib.modeling.backbone import Backbone
__all__ = ["DynamicNet4Seg", "SemSegDecoderHead", "BudgetConstraint"]

class DynamicNet4Seg(nn.Module):
    """
    This module implements Dynamic Network for Semantic Segmentation.
    """
    def __init__(self, config):
        super().__init__()
        self.constrain_on = config['MODEL']['BUDGET']['CONSTRAIN']
        self.unupdate_rate = config['MODEL']['BUDGET']['UNUPDATE_RATE']
        self.device = torch.device(config['MODEL']['DEVICE'])
        self.backbone = build_backbone(config)
        self.sem_seg_head = build_sem_seg_head(
            config, self.backbone.output_shape())
        pixel_mean = torch.Tensor(config['MODEL']['PIXEL_MEAN']).to(self.device).view(
            -1, 1, 1)
        pixel_std = torch.Tensor(config['MODEL']['PIXEL_STD']).to(self.device).view(
            -1, 1, 1)
        self.normalizer = lambda x: (x - pixel_mean) / pixel_std
        self.budget_constrint = BudgetConstraint(config)
        self.to(self.device)

    def forward(self, images, targets=None, step_rate=0.0):
        """
        Args:
            batched_inputs: a list, batched outputs of :class:`DatasetMapper` .
                Each item in the list contains the inputs for one image.
            step_rate: a float, calculated by current_step/total_step,
                This parameter is used for Scheduled Drop Path.
        For now, each item in the list is a dict that contains:
            image: Tensor, image in (C, H, W) format.
            sem_seg: semantic segmentation ground truth
            Other information that's included in the original dicts, such as:
                "height", "width" (int): the output resolution of the model, used in inference.
                    See :meth:`postprocess` for details.
        Returns:
            list[dict]: Each dict is the output for one input image.
                The dict contains one key "sem_seg" whose value is a
                Tensor of the output resolution that represents the
                per-pixel segmentation prediction.
        """
        features, expt_flops, real_flops, gate_list, gate_weight_list= self.backbone(images, step_rate)

        results, losses = self.sem_seg_head(features, targets)
        # calculate flops
        real_flops += self.sem_seg_head.flops
        flops = {'real_flops': real_flops, 'expt_flops': expt_flops}
        # use budget constraint for training
        if self.training and targets is not None:
            if self.constrain_on and step_rate >= self.unupdate_rate:
                warm_up_rate = min(
                    1.0, (step_rate - self.unupdate_rate) / 0.02
                )
                loss_budget = self.budget_constrint(
                    expt_flops, warm_up_rate=warm_up_rate
                )
                losses.update({'loss_budget': loss_budget})
            return losses, flops
        return results

class SemSegDecoderHead(nn.Module):
    """
    This module implements simple decoder head for Semantic Segmentation.
    It creats decoder on top of the dynamic backbone.
    """
    def __init__(self, config, input_shape: Dict[str, ShapeSpec]):
        super().__init__()
        # fmt: off
        self.in_features = config['MODEL']['SEM_SEG_HEAD']['IN_FEATURES']
        feature_strides = {k: v.stride for k, v in input_shape.items()}  # noqa:F841
        feature_channels = {k: v.channels for k, v in input_shape.items()}
        feature_resolution = {
            k: np.array([v.height, v.width])
            for k, v in input_shape.items()
        }
        self.ignore_value = config['MODEL']['SEM_SEG_HEAD']['IGNORE_VALUE']
        num_classes = config['MODEL']['SEM_SEG_HEAD']['NUM_CLASSES']
        norm = config['MODEL']['SEM_SEG_HEAD']['NORM']
        self.loss_weight = config['MODEL']['SEM_SEG_HEAD']['LOSS_WEIGHT']
        self.cal_flops = config['MODEL']['CAL_FLOPS']
        self.real_flops = 0.0
        # fmt: on

        self.layer_decoder_list = nn.ModuleList()
        # set affine in BatchNorm
        if 'Sync' in norm:
            affine = True
        else:
            affine = False
        # use simple decoder
        for _feat in self.in_features:
            res_size = feature_resolution[_feat]
            in_channel = feature_channels[_feat]
            if _feat == 'layer_0':
                out_channel = in_channel
            else:
                out_channel = in_channel // 2
            conv_1x1 = Conv2d(in_channel,
                              out_channel,
                              kernel_size=1,
                              stride=1,
                              padding=0,
                              bias=False,
                              norm=get_norm(norm, out_channel),
                              activation=nn.ReLU())
            self.real_flops += cal_op_flops.count_ConvBNReLU_flop(
                res_size[0],
                res_size[1],
                in_channel,
                out_channel, [1, 1],
                is_affine=affine)
            self.layer_decoder_list.append(conv_1x1)
        # using Kaiming init
        for layer in self.layer_decoder_list:
            weight_init.kaiming_init_module(layer, mode='fan_in')
        in_channel = feature_channels['layer_0']
        # the output layer
        self.predictor = Conv2d(in_channels=in_channel,
                                out_channels=num_classes,
                                kernel_size=3,
                                stride=1,
                                padding=1)
        self.real_flops += cal_op_flops.count_Conv_flop(
            feature_resolution['layer_0'][0], feature_resolution['layer_0'][1],
            in_channel, num_classes, [3, 3])
        # using Kaiming init
        weight_init.kaiming_init_module(self.predictor, mode='fan_in')

    def forward(self, features, targets=None):
        pred, pred_output = None, None
        for _index in range(len(self.in_features)):
            out_index = len(self.in_features) - _index - 1
            out_feat = features[self.in_features[out_index]]

            if isinstance(out_feat, float):
                continue

            if out_index <= 2 and pred is not None:
                out_feat = pred + out_feat

            pred = self.layer_decoder_list[out_index](out_feat)
            if out_index > 0:
                pred = F.interpolate(input=pred,
                                     scale_factor=2,
                                     mode='bilinear',
                                     align_corners=False)
            else:
                pred_output = pred
        # pred output
        pred_output = self.predictor(pred_output)
        pred_output = F.interpolate(input=pred_output,
                                    scale_factor=4,
                                    mode='bilinear',
                                    align_corners=False)

        if self.training and targets is not None:
            losses = {}
            losses["loss_sem_seg"] = (
                F.cross_entropy(
                    pred_output, targets, reduction="mean",
                    ignore_index=self.ignore_value
                ) * self.loss_weight
            )
            return [], losses
        else:
            return pred_output, {}

    @property
    def flops(self):
        return self.real_flops

class BudgetConstraint(nn.Module):
    """
    Given budget constraint to reduce expected inference FLOPs in the Dynamic Network.
    """
    def __init__(self, config):
        super().__init__()
        # fmt: off
        self.loss_weight = config['MODEL']['BUDGET']['LOSS_WEIGHT']
        self.loss_mu = config['MODEL']['BUDGET']['LOSS_MU']
        self.flops_all = config['MODEL']['BUDGET']['FLOPS_ALL']
        self.warm_up = config['MODEL']['BUDGET']['WARM_UP']
        # fmt: on

    def forward(self, flops_expt, warm_up_rate=1.0):
        if self.warm_up:
            warm_up_rate = min(1.0, warm_up_rate)
        else:
            warm_up_rate = 1.0
        losses = self.loss_weight * warm_up_rate * (
            (flops_expt / self.flops_all - self.loss_mu)**2
        )
        return losses

def build_backbone(config, input_shape=None):
    """
    Build a backbone from `config['MODEL.BACKBONE.NAME`.

    Returns:
        an instance of :class:`Backbone`
    """
    if input_shape is None:
        input_shape = ShapeSpec(channels=len(config['MODEL']['PIXEL_MEAN']),
                                height=config['INPUT']['FIX_SIZE_FOR_FLOPS'][0],
                                width=config['INPUT']['FIX_SIZE_FOR_FLOPS'][1])

    backbone = build_dynamic_backbone(config, input_shape)
    assert isinstance(backbone, Backbone)
    return backbone

def build_sem_seg_head(config, input_shape=None):
    return SemSegDecoderHead(config, input_shape)


text:
from abstractliftoveriter import AbstractLiftoverIter
import csv, math, h5py, os
import numpy as np

class AbstractMetaIter(AbstractLiftoverIter):
    def __init__(self, name, nameLong, directory, matrices, metadata, qtl, brainRegionFilter):
        super().__init__(name, directory)
        self.matrices = matrices
        self.metadata = metadata
        self.qtl = qtl
        self.brainRegionFilter = brainRegionFilter
        self.keys = []
        self.nameLong = nameLong

    def writeHDF5Metadata(self, root, rows):
        if not len(self.keys):
            raise "Must read circular rnas before selecting metadata"

        keyToIndexFiltered = {}
        counter = 0
        for row in rows:
            index = row.getMeta(self.id)
            if index != -1:
                keyToIndexFiltered[self.keys[index]] = counter
                counter += 1

        experimentGroup = root.create_group(self.name)
        experimentGroup.attrs.create("name", self.nameLong)

        sampleToIndex = None
        sampleHeadings = None
        if self.metadata:
            samples = experimentGroup.create_group("samples")
            sampleHeadings = self._writeHDF5Columns(self.metadata, samples)
            sampleToIndex = {}
            for i in range(len(sampleHeadings)):
                sampleToIndex[sampleHeadings[i]] = i

        if len(self.matrices):
            matrixGroup = None
            if self.metadata:
                matrixGroup = experimentGroup.create_group("matrices")
                matrixGroup.attrs.create("order", [m["type"] for m in self.matrices])
            for i in range(len(self.matrices)):
                heading, mdata2 = self._getAsMatrix(os.path.join(self.directory, self.matrices[i]["path"]), keyToIndexFiltered, sampleToIndex)
                arr = np.array(mdata2, dtype="f4")

                if sampleHeadings: heading = sampleHeadings

                if matrixGroup:
                    if not "sample_id" in experimentGroup:
                        experimentGroup.attrs.create("sample_id", np.array([h.encode() for h in heading], dtype="S" + str(len(max(heading, key=len)))))
                    matrixGroup.create_dataset(self.matrices[i]["type"], data=arr, chunks=(min(arr.shape[0], math.floor(10000/len(heading))), arr.shape[1]), compression="gzip", compression_opts=9)

                if i == 0:
                    self._writeHDF5Scaled(arr, experimentGroup)

        if self.qtl:
            self._writeHDF5CSVTable(self.qtl, experimentGroup, keyToIndexFiltered, rows)

        if self.brainRegionFilter:
            experimentGroup.attrs.create("brainRegionFilter", self.brainRegionFilter)

    def _getAsMatrix(self, fileName, keyToIndexFiltered, sampleToIndex, noneType="NA"):
        #Writes matrix in overall row order
        heading = None
        lines = [None] * len(keyToIndexFiltered)
        for line in csv.reader(open(fileName, 'r'), delimiter=','):
            if not heading: heading = line[1:]
            else:
                index = keyToIndexFiltered.get(line[0], -1)
                if index >= 0:
                    fromLine = line[1:]
                    if sampleToIndex:
                        toLine = [-1] * len(sampleToIndex)
                        for i in range(len(fromLine)):
                            subIndex = sampleToIndex.get(heading[i], -1)
                            if subIndex != -1:
                                toLine[subIndex] = fromLine[i]
                        lines[index] = toLine
                    else:
                        lines[index] = fromLine

        mdata2 = []
        for line in lines:
            mdata2.append(tuple([(line[i] if line[i] != noneType else -1.0) for i in range(len(line))]))
        return heading, mdata2

    #Only call this on cpm
    def _writeHDF5Scaled(self, arr, idGroup):
        colMeans = [np.mean([arr[i][j] for i in range(len(arr))]) for j in range(len(arr[0]))]
        colSTDs = [np.std([arr[i][j] for i in range(len(arr))]) for j in range(len(arr[0]))]

        logs = [np.mean([math.log2(abs(arr[i][j]) + 0.01) for j in range(len(arr[0]))]) for i in range(len(arr))]
        logMean = np.mean(logs)
        logSd = np.std(logs)
        arr2 = np.array([((x - logMean) / logSd) for x in logs], dtype="f4")
        idGroup.create_dataset("scaled", data=arr2, compression="gzip", compression_opts=9)

    def _writeHDF5Columns(self, fileName, hdf5Group, noneType="NA"):
        heading = []
        isHeading = True
        lines = []
        for line in csv.reader(open(os.path.join(self.directory, fileName), 'r'), delimiter=','):
            if isHeading:
                heading = line
                isHeading = False
            else:
                lines.append(line)
        heading[0] = "circ_id"

        allTypes = [int, float, str]
        allDefaults = [0, 0.0, ""]
        allTypesNp = ["i4", "f4", "S"]
        for i in range(1, len(lines[0])):
            for k in range(len(allTypes)):
                try:
                    # Attempt to parse all as this type
                    colType = allTypes[k]
                    values = [colType(lines[j][i]) for j in range(
                        len(lines)) if lines[j][i] != noneType]

                    # No exception so correct type, fix values for hdf5
                    colTypeNp = allTypesNp[k] + (str(len(max(values, key=len)))if allTypes[k] == str else "")
                    for m in range(len(values)):
                        if m == noneType: values[m] = allDefaults[k]
                        elif colType == str: values[m] = values[m].encode()

                    # Write dataset
                    arr = np.array(values, dtype=colTypeNp)
                    hdf5Group.create_dataset(heading[i], data=arr, chunks=arr.shape, compression="gzip", compression_opts=9)
                    break
                except:
                    continue
            else:
                raise("ERROR no type resolved for " + heading[i])
        hdf5Group.attrs.create("order", [heading[i] for i in range(1, len(heading))])

        return [lines[i][0] for i in range(len(lines))] #Return headings to copy this order

    def _writeHDF5CSVTable(self, fileName, hdf5Group, keyToIndexFiltered, rows, noneType="NA"):
        heading = []
        lines = [[] for i in range(len(keyToIndexFiltered))]

        qtlIter = open(os.path.join(self.directory, fileName), 'r').readlines().__iter__()
        heading = next(qtlIter).rstrip("\n").replace('"', '').split(",", 1)[1]

        #Now we know where to put everything that comes in
        for qtl in qtlIter:
            split = qtl.rstrip("\n").replace('"', '').split(',', 1)
            index = keyToIndexFiltered.get(split[0], -1)
            if index >= 0:
                lines[index].append(split[1])
        fixed = 0
        j = 0
        qtlIndices = [None] * len(lines)
        for i in range(len(lines)):
            qtlIndices[i] = fixed if len(lines[i]) else -1
            fixed += len(lines[i])

        hdf5Group.create_dataset("indices", data=np.array(qtlIndices, dtype="i4"), compression="gzip", compression_opts=9)

        mdata = []
        for i in range(len(lines)):
            if len(lines[i]) > 0:
                for subline in lines[i]:
                    mdata.append(subline)

        arr = np.array(mdata, dtype="S" + str(len(max(mdata, key=len))))
        ds = hdf5Group.create_dataset("QTL", data=arr, compression="gzip", compression_opts=9)
        ds.attrs.create("heading", heading)

max_stars_repo_name: Voineagulab/NeuroCirc

max_issues_repo_name: Voineagulab/NeuroCirc

max_forks_repo_name: Voineagulab/NeuroCirc


text:
import numpy as np
import matplotlib.pyplot as plt
import keras
from keras.datasets import mnist
from keras.models import Sequential
from keras.layers import Dense
from keras.optimizers import Adam
from keras.utils.np_utils import to_categorical
from keras.layers import Dropout, Flatten
from keras.layers.convolutional import Conv2D, MaxPooling2D
import random
import pickle
import pandas as pd

from keras.callbacks import LearningRateScheduler, ModelCheckpoint

np.random.seed(0)
with open('german-traffic-signs/train.p', 'rb') as f:
    train_data = pickle.load(f)
with open('german-traffic-signs/valid.p', 'rb') as f:
    val_data = pickle.load(f)
with open('german-traffic-signs/test.p', 'rb') as f:
    test_data = pickle.load(f)

print(type(train_data))
#print(train_data)

X_train, y_train = train_data['features'], train_data['labels']
X_val, y_val = val_data['features'], val_data['labels']
X_test, y_test = test_data['features'], test_data['labels']

print(X_train.shape)
print(X_val.shape)
print(X_test.shape)

assert(X_train.shape[0]==y_train.shape[0])    #The number of images is not equal to the number of labels
assert(X_val.shape[0]==y_val.shape[0])     #The number of images is not equal to the number of labels
assert(X_test.shape[0]==y_test.shape[0])     #The number of images is not equal to the number of labels
assert(X_train.shape[1:]==(32, 32, 3))    # The dimensions of the images are not 32x32x3
assert(X_val.shape[1:]==(32, 32, 3))      # The dimensions of the images are not 32x32x3
assert(X_test.shape[1:]==(32, 32, 3))   # The dimensions of the images are not 32x32x3

data = pd.read_csv('german-traffic-signs/signnames.csv')
print(data)
num_of_samples = []
cols = 5
num_classes = 43

fig, axs = plt.subplots(nrows=num_classes, ncols=cols, figsize=(5, 50))
fig.tight_layout()
#plt.show()

for i in range(cols):
    for j, row in data.iterrows():          # (index, Series)
        x_selected = X_train[y_train==j]
        axs[j][i].imshow(x_selected[random.randint(0, (len(x_selected) -1)), :, :], cmap=plt.get_cmap('gray'))
        axs[j][i].axis('off')
        if i == 2:
            axs[j][i].set_title(str(j) + '-' + row['SignName'])
            num_of_samples.append(len(x_selected))
plt.show()

print(num_of_samples)
plt.figure(figsize=(12, 4))
plt.bar(range(0, num_classes), num_of_samples)
plt.title('Distribution of the training dataset')
plt.xlabel('Class number')
plt.ylabel('Number of images')
plt.show()

import cv2

plt.imshow(X_train[1000])
plt.axis('off')
print(X_train[1000].shape)
print(y_train[1000])

def grayscale(img):
    img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    return img

img = grayscale(X_train[1000])
plt.imshow(img)
plt.axis('off')
print(img.shape)

def equalize(img):
    img = cv2.equalizeHist(img)
    return img

img = equalize(img)

plt.imshow(img)
plt.axis('off')
print(img.shape)

def preprocess(img):
    img = grayscale(img)
    img = equalize(img)
    img = img/255
    return img

X_train = np.array(list(map(preprocess, X_train)))
X_test = np.array(list(map(preprocess, X_test)))
X_val = np.array(list(map(preprocess, X_val)))

plt.imshow(X_train[random.randint(0, len(X_train)-1)])
plt.axis('off')
print(X_train.shape)

X_train = X_train.reshape(34799, 32, 32, 1)
X_test = X_test.reshape(12630, 32, 32, 1)
X_val = X_val.reshape(4410, 32, 32, 1)

from keras.preprocessing.image import ImageDataGenerator

datagen = ImageDataGenerator(width_shift_range=0.1,
                             height_shift_range=0.1,
                             zoom_range=0.2,
                             shear_range=0.1,
                             rotation_range=10.)

datagen.fit(X_train)

batches = datagen.flow(X_train, y_train, batch_size=15)
X_batch, y_batch = next(batches)

fig, axs = plt.subplots(1, 15, figsize=(20, 5))
fig.tight_layout()

for i in range(15):
    axs[i].imshow(X_batch[i].reshape(32, 32))
    axs[i].axis('off')

print(X_batch.shape)

y_train = to_categorical(y_train, 43)
y_test = to_categorical(y_test, 43)
y_val = to_categorical(y_val, 43)

def modified_model():
    model = Sequential()
    model.add(Conv2D(60, (5, 5), input_shape=(32, 32, 1), activation='relu'))
    model.add(Conv2D(60, (5, 5), activation='relu'))
    model.add(MaxPooling2D(pool_size=(2, 2)))

    model.add(Conv2D(30, (3, 3), activation='relu'))
    model.add(Conv2D(30, (3, 3), activation='relu'))
    model.add(MaxPooling2D(pool_size=(2, 2)))

    model.add(Flatten())
    model.add(Dense(500, activation='relu'))
    model.add(Dropout(0.5))
    model.add(Dense(43, activation='softmax'))

    model.compile(Adam(lr=0.001), loss='categorical_crossentropy', metrics=['accuracy'])
    return model

model = modified_model()

print(model.summary())

history = model.fit_generator(datagen.flow(X_train, y_train, batch_size=50),
                              steps_per_epoch=2000,
                              epochs=10,
                              validation_data=(X_val, y_val),
                              shuffle=1)

plt.plot(history.history['loss'])
plt.plot(history.history['val_loss'])
plt.legend(['loss', 'val_loss'])
plt.title('Loss')
plt.xlabel('epoch')

plt.plot(history.history['acc'])
plt.plot(history.history['val_acc'])
plt.legend(['acc', 'val_acc'])
plt.title('Accuracy')
plt.xlabel('epoch')

score = model.evaluate(X_test, y_test, verbose=0)
print(type(score))
print('Test score:', score[0])
print('Test accuracy:', score[1])

import requests
from PIL import Image

url = 'https://c8.alamy.com/comp/G667W0/road-sign-speed-limit-30-kmh-zone-passau-bavaria-germany-G667W0.jpg'
response = requests.get(url, stream=True)
print(response)
img = Image.open(response.raw)
plt.imshow(img, cmap=plt.get_cmap('gray'))

img = np.asarray(img)
print(img.shape)
img = cv2.resize(img, (32, 32))
img = preprocess(img)
plt.imshow(img, cmap=plt.get_cmap('gray'))
print(img.shape)

img = img.reshape(1, 32, 32, 1)
print(img.shape)

prediction = model.predict_classes(img)
print('prediction digit:', str(prediction))

max_stars_repo_path: Traffic Signs/Road Symbols.py

max_issues_repo_path: Traffic Signs/Road Symbols.py

max_forks_repo_path: Traffic Signs/Road Symbols.py


text:
from retinanet.dataset import Ring_Cell_all_dataset
from tqdm import tqdm
import numpy as np
import torch
# from lib.nms.pth_nms import pth_nms
from lib_new.nms.gpu_nms import gpu_nms
from metric import compute_overlap
import matplotlib.pyplot as plt
import matplotlib
import os
import glob

result_str = ''

for test_fold in range(0,4):
	test_dataset = Ring_Cell_all_dataset('/data/sqy/code/miccai2019/train_test_4/test_{}.txt'.format(test_fold))

	pred_dict_box = {}
	pred_dict_score = {}
	pred_dict_max_valid_recall = {}

	csv_path = './test_result_new/retinanet_resnet18_round0_fold_{}_weight_loss_1_on_test_data_latest.csv'.format(test_fold)

	with open(csv_path, 'r') as f:
		lines = f.readlines()
		for line in lines:
			line = line[:-1]
			line = line.split(',')
			image_name = line[0]
			pred_dict_box[image_name] = []
			pred_dict_score[image_name] = []
			if len(line[1]) != 0:
				preds = line[1].split(';')[:-1]
				for pred in preds:
					pred = pred.split(' ')
					box = []
					for elemet in pred[:-1]:
						box.append(float(elemet))
					pred_dict_box[image_name].append(box)
					pred_dict_score[image_name].append(float(pred[-1]))

	result_dict = {}
	nms_threshold = 0.4
	scores = np.zeros((0,))
	iou_threshold = 0.3
	scores_all_prediction = []

	for i, (image, bbox, image_, image_name) in enumerate(tqdm(test_dataset)):
		result_dict[image_name] = []
		detected_annotations = []
		scores_all_prediction.extend(pred_dict_score[image_name])

		false_positives = np.zeros((0,))
		true_positives = np.zeros((0,))
		scores = np.zeros((0,))
		num_annotations = len(bbox)

		if len(bbox) != 0:
			# pos image
			pred_scores = pred_dict_score[image_name]
			detections = pred_dict_box[image_name]
			annotations = np.array(bbox)
			for j, d in enumerate(detections):
				score = pred_scores[j]
				scores = np.append(scores, score)
				overlaps = compute_overlap(np.expand_dims(d, axis=0), annotations)
				assigned_annotation = np.argmax(overlaps, axis=1)
				max_overlap = overlaps[0, assigned_annotation]

				if max_overlap >= iou_threshold and assigned_annotation not in detected_annotations:
					detected_annotations.append(assigned_annotation)
					gt_box = bbox[assigned_annotation[0]]
					gt_box.append(score)
					result_dict[image_name].append(gt_box)
					false_positives = np.append(false_positives, 0)
					true_positives = np.append(true_positives, 1)
				else:
					false_positives = np.append(false_positives, 1)
					true_positives = np.append(true_positives, 0)

			indices = np.argsort(-scores)
			scores = scores[indices]
			false_positives = false_positives[indices]
			true_positives = true_positives[indices]

			false_positives = np.cumsum(false_positives)
			true_positives = np.cumsum(true_positives)

			recall = true_positives / num_annotations
			precision = true_positives / np.maximum(true_positives + false_positives, np.finfo(np.float64).eps)

			index_record = np.where(precision >= 0.2)[0]

			if index_record.shape[0] != 0:
				index_record = index_record[-1]
				# recall, precision, FPs when precision is 0.2
				recall_record = recall[:index_record + 1]
				precision_record = precision[:index_record + 1]
				score_record = scores[:index_record + 1]
			else:
				recall_record = [0]
				precision_record = [0]
				score_record = [0]

			pred_dict_max_valid_recall[image_name] = recall_record[-1]

		for j in range(len(bbox)):
			if j not in detected_annotations:
				gt_box = bbox[j]
				gt_box.append(0)
				result_dict[image_name].append(gt_box)

	for image_name in result_dict:
		result_str += image_name
		result_str += ','
		results = result_dict[image_name]
		for result in results:
			box = result
			for element in box:
				result_str += str(element)
				result_str += ' '
			result_str += ';'

		if image_name not in pred_dict_max_valid_recall:
			result_str += '-1'
		else:
			result_str += str(pred_dict_max_valid_recall[image_name])

		result_str += '\n'

result_csv = './bbox/retinanet_resnet18_training_data_with_confidence_score_of_whole_image.csv'

with open(result_csv, 'w') as f:
	f.write(result_str)


text:
#!/usr/bin/env python

# import modules
import rospy
import tf
from kuka_arm.srv import *
from trajectory_msgs.msg import JointTrajectory, JointTrajectoryPoint
from geometry_msgs.msg import Pose
from mpmath import *
from sympy import *

def handle_calculate_IK(req):
    rospy.loginfo("Received %s eef-poses from the plan" % len(req.poses))
    if len(req.poses) < 1:
        print "No valid poses received"
        return -1
    else:

        ### Your FK code here
        # Create symbols
	d1, d2, d3, d4, d5, d6, d7 = symbols('d1:8')
        a0, a1, a2, a3, a4, a5, a6 = symbols('a0:7')
   	alpha0, alpha1, alpha2, alpha3, alpha4, alpha5, alpha6 = symbols('alpha0:7')
    	q1, q2, q3, q4, q5, q6, q7 = symbols('q1:8')
	#
	# Create Modified DH parameters
	s = {alpha0:       0, a0:      0, d1:  0.75, q1:      q1,
    	    alpha1:    -pi/2, a1:   0.35, d2:     0, q2: q2-pi/2,
     	    alpha2:        0, a2:   1.25, d3:     0, q3:	  q3,
     	    alpha3:    -pi/2, a3: -0.054, d4:   1.5, q4:	  q1,
     	    alpha4:     pi/2, a4:      0, d5:     0, q5:	  q5,
     	    alpha5:    -pi/2, a5:      0, d6:     0, q6:	  q6,
     	    alpha6:        0, a6:      0, d7: 0.303, q7:      0}

	#
	# Define Modified DH Transformation matrix
	def Rot_x(q):
            Rx = Matrix([[ 1,      0,       0],
                      [ 0, cos(q), -sin(q)],
                      [ 0, sin(q),  cos(q)]])
            return Rx

    	def Rot_y(q):
            Ry = Matrix([[  cos(q), 0, sin(q)],
                      [       0, 1,      0],
                      [ -sin(q), 0, cos(q)]])
            return Ry

    	def Rot_z(q):
            Rz = Matrix([[ cos(q), -sin(q), 0],
                      [ sin(q),  cos(q), 0],
                      [      0,       0, 1]])
            return Rz

    	def TF_Matrix(alpha, a, d, q):
	    TF = Matrix([[             cos(q),            -sin(q),            0,              a],
	 	[ sin(q)*cos(alpha), cos(q)*cos(alpha), -sin(alpha), -sin(alpha)*d],
		[ sin(q)*sin(alpha), cos(q)*sin(alpha),  cos(alpha),  cos(alpha)*d],
		[                   0,                   0,            0,               1]])
    	    return TF
	# Conversion factors between radians and degrees
    	rtd = 180 / pi
    	dtr = pi / 180
	#
	# Create individual transformation matrices
	 # Create individual transformation matrices, for now
    	T0_1=TF_Matrix(alpha0, a0, d1, q1).subs(s)
    	T1_2=TF_Matrix(alpha1, a1, d2, q2).subs(s)
    	T2_3=TF_Matrix(alpha2, a2, d3, q3).subs(s)
    	T3_4=TF_Matrix(alpha3, a3, d4, q4).subs(s)
    	T4_5=TF_Matrix(alpha4, a4, d5, q5).subs(s)
    	T5_6=TF_Matrix(alpha5, a5, d6, q6).subs(s)
    	T6_G=TF_Matrix(alpha6, a6, d7, q7).subs(s)

    	T0_G=T0_1 * T1_2 * T2_3 * T3_4 * T4_5 * T5_6 * T6_G
	R0_3 = T0_1[0:3, 0:3]*T1_2[0:3, 0:3]*T2_3[0:3, 0:3]
	# Extract rotation matrices from the transformation matrices
	#
	#
        ###

        # Initialize service response
        joint_trajectory_list = []
        for x in xrange(0, len(req.poses)):
            # IK code starts here
            joint_trajectory_point = JointTrajectoryPoint()

	    # Extract end-effector position and orientation from request
	    # px,py,pz = end-effector position
	    # roll, pitch, yaw = end-effector orientation
            px = req.poses[x].position.x
            py = req.poses[x].position.y
            pz = req.poses[x].position.z

            (roll, pitch, yaw) = tf.transformations.euler_from_quaternion(
                [req.poses[x].orientation.x, req.poses[x].orientation.y,
                    req.poses[x].orientation.z, req.poses[x].orientation.w])

            ### Your IK code here
	    # Compensate for rotation discrepancy between DH parameters and Gazebo
	    r, p, y = symbols('r p y')
   	    #Solution is in DH Reference freme, to convert it to URDF file  reference frame we rortate 180 deg around Z and -90 around Y
	    Rot_corr=Rot_z(180*dtr)*Rot_y(-90*dtr)
    	    R_EE = Rot_z(yaw)*Rot_y(pitch)*Rot_x(roll)*Rot_corr

   	    EE=Matrix([[px],[py],[pz]]) ##EE point location

    	    ##WC is only a translation difference from J EE_to J6 along, the Z axis of the local reference frame
    	    #wrist center
    	    WC= EE - (0.303)*R_EE[:,2]
	    #
	    # Calculate joint angles using Geometric IK method
	    theta1=atan2(WC[1],WC[0])
    	    side_a=1.50
    	    side_b_xy=sqrt(WC[0]*WC[0]+WC[1]*WC[1])-0.35
    	    side_b_z=WC[2]-0.75
   	    side_b=sqrt(pow((side_b_xy),2) + pow((side_b_z), 2))
   	    side_c=1.25
    	    #Angles using cosine laws
    	    angle_a=acos((   side_b * side_b + side_c * side_c - side_a * side_a )/( 2 * side_b * side_c ))
    	    angle_b=acos(( - side_b * side_b + side_c * side_c + side_a * side_a )/( 2 * side_a * side_c ))
    	    angle_c=acos((   side_b * side_b - side_c * side_c + side_a * side_a )/( 2 * side_a * side_b ))
    	    ### Thetas
    	    theta2= pi/2 - angle_a - atan2(side_b_z, side_b_xy)
    	    theta3= pi/2 - angle_b + 0.036 #sag in link 4
    	    #multiplying the the R in T trans     matrieces. (picking the first 3 rows and 3 columns)
    	    R0_3 = R0_3.evalf(subs={q1: theta1, q2: theta2, q3: theta3})
    	    R3_6 = R0_3.T * R_EE

    	    theta4 = atan2(R3_6[2, 2], -R3_6[0, 2])
    	    theta5 = atan2(sqrt(R3_6[0, 2]**2+R3_6[2, 2]**2), R3_6[1, 2])
    	    theta6 = atan2(-R3_6[1, 1], R3_6[1, 0])

    	    # In the next line replace theta1,theta2...,theta6 by your joint angle variables
	    #
            ###

            # Populate response for the IK request
            # In the next line replace theta1,theta2...,theta6 by your joint angle variables
	    joint_trajectory_point.positions = [theta1, theta2, theta3, theta4, theta5, theta6]
	    joint_trajectory_list.append(joint_trajectory_point)

        rospy.loginfo("length of Joint Trajectory List: %s" % len(joint_trajectory_list))
        return CalculateIKResponse(joint_trajectory_list)

def IK_server():
    # initialize node and declare calculate_ik service
    rospy.init_node('IK_server')
    s = rospy.Service('calculate_ik', CalculateIK, handle_calculate_IK)
    print "Ready to receive an IK request"
    rospy.spin()

if __name__ == "__main__":
    IK_server()

max_stars_repo_name: githubpramod/Robotic-Arm-Pick-Place

max_issues_repo_name: githubpramod/Robotic-Arm-Pick-Place

max_forks_repo_name: githubpramod/Robotic-Arm-Pick-Place

reason: from sympy,from mpmath


text:
"""GaussianKDE module."""

import numpy as np
from scipy.special import ndtr
from scipy.stats import gaussian_kde

from copulas import EPSILON, random_state, store_args, validate_random_state
from copulas.optimize import bisect, chandrupatla
from copulas.univariate.base import BoundedType, ParametricType, ScipyModel

class GaussianKDE(ScipyModel):
    """A wrapper for gaussian Kernel density estimation.

    It was implemented in scipy.stats toolbox. gaussian_kde is slower than statsmodels
    but allows more flexibility.

    When a sample_size is provided the fit method will sample the
    data, and mask the real information. Also, ensure the number of
    entries will be always the value of sample_size.

    Args:
        sample_size(int): amount of parameters to sample
    """

    PARAMETRIC = ParametricType.NON_PARAMETRIC
    BOUNDED = BoundedType.UNBOUNDED
    MODEL_CLASS = gaussian_kde

    @store_args
    def __init__(self, sample_size=None, random_state=None, bw_method=None, weights=None):
        self.random_state = validate_random_state(random_state)
        self._sample_size = sample_size
        self.bw_method = bw_method
        self.weights = weights

    def _get_model(self):
        dataset = self._params['dataset']
        self._sample_size = self._sample_size or len(dataset)
        return gaussian_kde(dataset, bw_method=self.bw_method, weights=self.weights)

    def _get_bounds(self):
        X = self._params['dataset']
        lower = np.min(X) - (5 * np.std(X))
        upper = np.max(X) + (5 * np.std(X))

        return lower, upper

    def probability_density(self, X):
        """Compute the probability density for each point in X.

        Arguments:
            X (numpy.ndarray):
                Values for which the probability density will be computed.
                It must have shape (n, 1).

        Returns:
            numpy.ndarray:
                Probability density values for points in X.

        Raises:
            NotFittedError:
                if the model is not fitted.
        """
        self.check_fit()
        return self._model.evaluate(X)

    @random_state
    def sample(self, n_samples=1):
        """Sample values from this model.

        Argument:
            n_samples (int):
                Number of values to sample

        Returns:
            numpy.ndarray:
                Array of shape (n_samples, 1) with values randomly
                sampled from this model distribution.

        Raises:
            NotFittedError:
                if the model is not fitted.
        """
        self.check_fit()
        return self._model.resample(size=n_samples)[0]

    def cumulative_distribution(self, X):
        """Compute the cumulative distribution value for each point in X.

        Arguments:
            X (numpy.ndarray):
                Values for which the cumulative distribution will be computed.
                It must have shape (n, 1).

        Returns:
            numpy.ndarray:
                Cumulative distribution values for points in X.

        Raises:
            NotFittedError:
                if the model is not fitted.
        """
        self.check_fit()
        X = np.array(X)
        stdev = np.sqrt(self._model.covariance[0, 0])
        lower = ndtr((self._get_bounds()[0] - self._model.dataset) / stdev)[0]
        uppers = ndtr((X[:, None] - self._model.dataset) / stdev)
        return (uppers - lower).dot(self._model.weights)

    def percent_point(self, U, method='chandrupatla'):
        """Compute the inverse cumulative distribution value for each point in U.

        Arguments:
            U (numpy.ndarray):
                Values for which the cumulative distribution will be computed.
                It must have shape (n, 1) and values must be in [0,1].
            method (str):
                Whether to use the `chandrupatla` or `bisect` solver.

        Returns:
            numpy.ndarray:
                Inverse cumulative distribution values for points in U.

        Raises:
            NotFittedError:
                if the model is not fitted.
        """
        self.check_fit()

        if len(U.shape) > 1:
            raise ValueError(f'Expected 1d array, got {(U, )}.')

        if np.any(U > 1.0) or np.any(U < 0.0):
            raise ValueError('Expected values in range [0.0, 1.0].')

        is_one = U >= 1.0 - EPSILON
        is_zero = U <= EPSILON
        is_valid = ~(is_zero | is_one)

        lower, upper = self._get_bounds()

        def _f(X):
            return self.cumulative_distribution(X) - U[is_valid]

        X = np.zeros(U.shape)
        X[is_one] = float('inf')
        X[is_zero] = float('-inf')
        if is_valid.any():
            lower = np.full(U[is_valid].shape, lower)
            upper = np.full(U[is_valid].shape, upper)
            if method == 'bisect':
                X[is_valid] = bisect(_f, lower, upper)
            else:
                X[is_valid] = chandrupatla(_f, lower, upper)

        return X

    def _fit_constant(self, X):
        sample_size = self._sample_size or len(X)
        constant = np.unique(X)[0]
        self._params = {
            'dataset': [constant] * sample_size,
        }

    def _fit(self, X):
        if self._sample_size:
            X = gaussian_kde(X, bw_method=self.bw_method,
                             weights=self.weights).resample(self._sample_size)
        self._params = {
            'dataset': X.tolist()
        }
        self._model = self._get_model()

    def _is_constant(self):
        return len(np.unique(self._params['dataset'])) == 1

    def _extract_constant(self):
        return self._params['dataset'][0]

    def _set_params(self, params):
        """Set the parameters of this univariate.

        Args:
            params (dict):
                Parameters to recreate this instance.
        """
        self._params = params.copy()
        if self._is_constant():
            constant = self._extract_constant()
            self._set_constant_value(constant)
        else:
            self._model = self._get_model()

max_stars_repo_name: DAI-Lab/Copulas

max_issues_repo_name: DAI-Lab/Copulas

max_forks_repo_name: DAI-Lab/Copulas

reason: import numpy,from scipy


text:
# https://arxiv.org/pdf/1708.03080.pdf
# https://www.google.com/url?sa=t&rct=j&q=&esrc=s&source=web&cd=&cad=rja&uact=8&ved=2ahUKEwjFlNOipOL1AhXnjokEHU6LAnUQFnoECA0QAw&url=https%3A%2F%2Fwww.sciencedirect.com%2Fscience%2Farticle%2Fpii%2FS2352146520304804%2Fpdf%3Fmd5%3D7efbb0aed0d1a2cdec7d2fb6d040f632%26pid%3D1-s2.0-S2352146520304804-main.pdf&usg=AOvVaw3A5ZPzSHEk03XlNIj-TR3C

from math import *
from typing import List, Sequence
import matplotlib.pyplot as plt
from random import *
import numpy as np
from PIL import Image
import heapq
import cv2
import os
import json
import cProfile

AGENT_RADIUS = 0.2  # size of agent (m)
AGENT_MAX_ANGLE_CHANGE = pi/2  # maximum agent direction change
AGENT_DISTANCE_WEIGHT = 0.4  # lower = more likely to slow down
AGENT_COMFORTABLE_SPACING = 1 # additional distance the agent likes to keep to others (m)
AGENT_SPACING_SHIFT = 0.8 # agent shift away from others (per step) (m) (NOT as a function of step size)
AGENT_STEP = 0.8  # step size of the agent (m)
AGENT_SEATING_PENALTY = 5 # scalar for the agent's dislike for seats (1 implies no avoidance of seats)
AGENT_PATH_LOOKAHEAD = 2 # number of steps ahead of current path location the agent tries to move to
AGENT_STUCK_THRESHOLD = 6 # how quickly the agent decides it is stuck and changes to a more direct pathing
AGENT_CROWD_THRESHOLD = 30 # how many people have to be present for the agent to decide "this is a crowd" and decrease movement checks
AGENT_THICK_CROWD_THRESHOLD = 52 # how many people have to be present for the agent to ignore shift calculations
EXIT_WIDTH = 2 # how close the people need to be to the exit to count as finished (path steps)
WORLD_FILE = "Stadium-map.png"
WORLD_EXTENTS = Image.open(WORLD_FILE).size

COLLISION_DIST_SQUARED = (2 * AGENT_RADIUS) * (2 * AGENT_RADIUS)
AVOIDANCE_DIST_SQUARED = (2 * AGENT_RADIUS + AGENT_COMFORTABLE_SPACING) * (2 * AGENT_RADIUS + AGENT_COMFORTABLE_SPACING)

values = {}
arcs = []
for theta in range(0, 5):
    arcs.append((0.5, (theta-2) / 2 * AGENT_MAX_ANGLE_CHANGE))
    for a in range(1, 5):
        direction = (theta-2) / 2 * AGENT_MAX_ANGLE_CHANGE
        distance = a / 4
        value = AGENT_DISTANCE_WEIGHT * a + (1 - AGENT_DISTANCE_WEIGHT) * (1 - abs(direction) / AGENT_MAX_ANGLE_CHANGE)
        values[(distance, direction)] = value

moves = [x for x in values]
moves.sort(key=lambda x: values[x], reverse=True)

class Agent:
    x_pos = 0
    y_pos = 0
    desired_step_size = 1.0
    desired_location = (0, 0)
    collision_grid_pos = (0, 0)
    desired_direction = 0
    path = None
    finished = False
    stuck_factor = 0
    stay_put_timer = 0

    def __init__(self, pos, speed, goal, time):
        self.x_pos, self.y_pos = pos
        self.desired_step_size = speed
        self.desired_location = goal
        self.stay_put_timer = time
        self.compute_path()

    def compute_path(self):
        global world_collisions
        path = world_collisions.get_path((floor(self.x_pos), floor(self.y_pos)), (floor(self.desired_location[0]), floor(self.desired_location[1])))
        if path:
            if path == self.path:
                self.stuck_factor += 1
            else:
                self.path = path
                self.stuck_factor = 0

    def update(self):
        if self.stay_put_timer > 0:
            self.stay_put_timer -= AGENT_STEP/1.4
            return

        assert not self.finished

        global agent_collisions
        global world_collisions

        self.compute_path()

        if len(self.path) <= EXIT_WIDTH:
            agent_collisions.deregister_member(self)
            self.finished = True
            return

        lookahead = min(len(self.path)-1,AGENT_PATH_LOOKAHEAD)
        if self.stuck_factor >= AGENT_STUCK_THRESHOLD:
            lookahead = 1
        self.desired_direction = atan2(self.path[lookahead][1]+0.5 - self.y_pos, self.path[lookahead][0]+0.5 - self.x_pos)

        nearby_people = agent_collisions.get_all_collision_considerations(self.collision_grid_pos)

        if len(nearby_people) > AGENT_CROWD_THRESHOLD and self.stuck_factor >= 10*AGENT_STUCK_THRESHOLD:
            self.desired_location = world_collisions.random_goal()
            self.stuck_factor = 0
            self.update()

        shift = self.shift_avoid_people(nearby_people)
        original_length = self.desired_step_size
        original_direction = self.desired_direction
        self.desired_step_size = shift[0]
        self.desired_direction = shift[1]

        move = self.find_next_move(nearby_people)
        self.x_pos, self.y_pos = self.next_position(move)
        agent_collisions.check_for_updates(self)

        self.desired_step_size = original_length
        self.desired_direction = original_direction

    def shift_avoid_people(self, people):
        # if len(people) > AGENT_THICK_CROWD_THRESHOLD:
        #     return [self.desired_step_size, self.desired_direction]

        global agent_collisions
        people_pos_sum = [0, 0]
        n_people = 0
        for obstacle in people:
            if self != obstacle:
                dx = (self.x_pos - obstacle.x_pos)
                dy = (self.y_pos - obstacle.y_pos)
                if dx*dx + dy*dy < AVOIDANCE_DIST_SQUARED:
                    people_pos_sum[0] += obstacle.x_pos
                    people_pos_sum[1] += obstacle.y_pos
                    n_people += 1
        if n_people > 0:
            people_center = [people_pos_sum[0] / n_people, people_pos_sum[1] / n_people]
            avoidance_vector = [self.x_pos - people_center[0], self.y_pos - people_center[1]]
            avoidance_dir = atan2(avoidance_vector[1], avoidance_vector[0])
            shift = [
                cos(avoidance_dir) * AGENT_SPACING_SHIFT,
                sin(avoidance_dir) * AGENT_SPACING_SHIFT
                ]
            optimal_vector = [
                cos(self.desired_direction) * self.desired_step_size,
                sin(self.desired_direction) * self.desired_step_size
                ]
            shifted_vector = [
                optimal_vector[0] + shift[0],
                optimal_vector[1] + shift[1]
                ]
            shifted_angle = atan2(shifted_vector[1], shifted_vector[0])
            shifted_length = sqrt(shifted_vector[0]**2 + shifted_vector[1]**2)
            shifted_length = max(min(shifted_length, self.desired_step_size), 0)
            return [self.desired_step_size, shifted_angle]
        return [self.desired_step_size, self.desired_direction]

    def next_position(self, move):
        # according to the paper, sin and cos are flipped here, but that doesn't make sense, soooooooo...
        x = move[0] * self.desired_step_size * cos(self.desired_direction + move[1])
        y = move[0] * self.desired_step_size * sin(self.desired_direction + move[1])
        return (self.x_pos + x, self.y_pos + y)

    def find_next_move(self, people):
        global moves
        global arcs
        illegal_directions = set()
        if len(people) >= AGENT_CROWD_THRESHOLD:
            for arc in arcs:
                n = self.next_position(arc)
                if not self.is_legal_move(n, people):
                    illegal_directions.add(arc[1])

        best_selection = (0.0, 0.0)
        for move in moves:
            if move[1] in illegal_directions:
                continue
            n = self.next_position(move)
            if self.is_legal_move(n, people):
                best_selection = move
                break
        return best_selection

    def is_legal_move(self, n, people):
        global agent_collisions
        global world_collisions
        if not world_collisions.is_valid_location(n):
            return False
        for obstacle in people:
            if self != obstacle:
                dx = (n[0] - obstacle.x_pos)
                dy = (n[1] - obstacle.y_pos)
                if dx*dx + dy*dy < COLLISION_DIST_SQUARED:
                    return False
        return True

class AgentCollisionManager:
    square_side = 1
    squares = []

    def __init__(self, x, y):
        self.squares = [[[] for _ in range(ceil(y/self.square_side)+1)] for _ in range(ceil(x/self.square_side)+1)]
        # +1 because buffer (avoids having to out-of-bounds check, saving time in get_all_collision_considerations)

    def check_for_updates(self, member: Agent):
        coord = [floor(member.x_pos/self.square_side),floor(member.y_pos/self.square_side)]
        if coord[0] != member.collision_grid_pos[0] and coord[1] != member.collision_grid_pos[1]:
            self.squares[member.collision_grid_pos[0]][member.collision_grid_pos[1]].remove(member)
            self.squares[coord[0]][coord[1]].append(member)
            member.collision_grid_pos = coord

    def register_member(self, member: Agent):
        coord = [floor(member.x_pos/self.square_side),floor(member.y_pos/self.square_side)]
        self.squares[coord[0]][coord[1]].append(member)
        member.collision_grid_pos = coord

    def deregister_member(self, member: Agent):
        self.squares[member.collision_grid_pos[0]][member.collision_grid_pos[1]].remove(member)

    def point_in_square(self, pos, index):
        bounds = self.get_square_bounds(index)
        return pos[0] >= bounds[0] and pos[1] >= bounds[1] and pos[0] < bounds[2] and pos[1] < bounds[3]

    def get_square_bounds(self, index):
        min_x = index[0] * self.square_side
        min_y = index[1] * self.square_side
        max_x = (index[0] + 1) * self.square_side
        max_y = (index[1] + 1) * self.square_side
        return [min_x, min_y, max_x, max_y]

    def get_all_collision_considerations(self, index):
        objects = []
        for dx in range(-1, 2):
            for dy in range(-1, 2):
                objects += self.squares[index[0]+dx][index[1]+dy]
        return objects

class WorldManager:
    illegal_pixels = set()
    seating_pixels = set()
    goals = []
    ideal_goals = {}
    paths = {}

    def __init__(self, image):
        arr = np.asarray(image)
        for x in range(len(arr[0])):
            for y in range(len(arr)):
                px = arr[(WORLD_EXTENTS[1]-1)-y, x]
                if px[0] == 0 and px[1] == 0 and px[2] == 0:
                    self.illegal_pixels.add((float(x), float(y)))
                if px[0] == 255 and px[1] == 0 and px[2] == 0:
                    self.seating_pixels.add((float(x), float(y)))
                if px[0] == 0 and px[1] == 255 and px[2] == 0:
                    self.goals.append((x, y))
        self.compute_ideal_goals()

    def read_paths(self):
        if os.path.exists(WORLD_FILE[:-4] + "-pathdata.json"):
            f = open(WORLD_FILE[:-4] + "-pathdata.json", "r")
            self.paths = json.loads(f.read())
            f.close()
        else:
            self.precompute_paths()

    def precompute_paths(self):
        for goal in self.goals:
            print("computing paths to:", goal)
            self.compute_paths(goal)
        print("done computing paths")
        with open(WORLD_FILE[:-4] + "-pathdata.json", "w") as f:
            f.write(json.dumps(self.paths))

    def is_valid_location(self, point):
        x = floor(point[0])
        y = floor(point[1])
        return (x, y) not in self.illegal_pixels

    def is_seating_location(self, point):
        x = floor(point[0])
        y = floor(point[1])
        return (x, y) in self.seating_pixels

    def compute_ideal_goals(self):
        if os.path.exists(WORLD_FILE[:-4] + "-exitdata.json"):
            f = open(WORLD_FILE[:-4] + "-exitdata.json", "r")
            self.ideal_goals = json.loads(f.read())
            f.close()
        else:
            goal_preferences = [[] for _ in range(len(self.goals))]
            for i in range(len(self.goals)):
                for pos in self.seating_pixels:
                    goal_preferences[i].append([(pos[0]-self.goals[i][0])**2 + (pos[1]-self.goals[i][1])**2, pos])
                goal_preferences[i].sort(key=lambda x: x[0])

            assigned = set()
            goal_considered = [0 for _ in range(len(self.goals))]
            i = 0
            while len(assigned) < len(self.seating_pixels):
                while goal_preferences[i][goal_considered[i]][1] in assigned:
                    goal_considered[i] +=1
                assigned.add(goal_preferences[i][goal_considered[i]][1])
                self.ideal_goals[str(goal_preferences[i][goal_considered[i]][1])] = self.goals[i]
                i = (i+1)%len(self.goals)

            with open(WORLD_FILE[:-4] + "-exitdata.json", "w") as f:
                f.write(json.dumps(self.ideal_goals))

    def smart_closest_goal(self, pos):
        return self.ideal_goals[str((float(floor(pos[0])),float(floor(pos[1]))))]

    def closest_goal(self, pos):
        best = 0
        distance = 1e100
        for i,g in enumerate(self.goals):
            if (g[0]-pos[0])**2 + (g[1]-pos[1])**2 < distance:
                distance = (g[0]-pos[0])**2 + (g[1]-pos[1])**2
                best = i
        return self.goals[best]

    def random_goal(self):
        return self.goals[randint(0,len(self.goals)-1)]

    def assign_goal(self, mode, pos):
        if mode == "random":
            return self.random_goal()
        elif mode == "closest":
            return self.closest_goal(pos)
        elif mode == "smartclosest":
            return self.smart_closest_goal(pos)

    def get_path(self, pos, goal):
        if str((pos[0],pos[1],goal[0],goal[1])) in self.paths:
            path = [pos]
            for _ in range(AGENT_PATH_LOOKAHEAD):
                p = self.paths[str((path[-1][0],path[-1][1],goal[0],goal[1]))]
                if p == None:
                    break
                path.append(p)
            return path
        else:
            print(str((pos[0],pos[1],goal[0],goal[1])))

    def compute_paths(self, goal):
        # Djikstra's algorithm
        # Computes every path to this goal
        frontier = [DijkstraNode(goal, None)]
        heapq.heapify(frontier)
        frontier_set = set(frontier)
        visited = set()

        while len(frontier) != 0:
            next_node = heapq.heappop(frontier)
            visited.add(next_node)
            frontier_set.remove(next_node)
            for neighbor in self.neighbors(next_node.pos):
                if not self.is_valid_location(neighbor):
                    continue
                n = DijkstraNode(neighbor, next_node)
                if n in visited:
                    continue
                if n not in frontier_set:
                    frontier_set.add(n)
                    heapq.heappush(frontier, n)
                else:
                    for e in frontier_set:
                        if e == n:
                            if e > n:
                                e.set_parent(next_node)
                            break

        for node in visited:
            if node.parent != None:
                self.paths[str((node.pos[0], node.pos[1], goal[0], goal[1]))] = node.parent.pos
            else:
                self.paths[str((node.pos[0], node.pos[1], goal[0], goal[1]))] = None

    def neighbors(self, pos):
        return [
            (pos[0] - 1, pos[1] - 1),
            (pos[0], pos[1] - 1),
            (pos[0] + 1, pos[1] - 1),
            (pos[0] - 1, pos[1]),
            (pos[0] + 1, pos[1]),
            (pos[0] - 1, pos[1] + 1),
            (pos[0], pos[1] + 1),
            (pos[0] + 1, pos[1] + 1)
        ]

class DijkstraNode:
    pos = (0, 0)
    parent = None
    cost = 0

    def __init__(self, pos, parent):
        self.pos = pos
        self.parent = parent
        if self.parent == None:
            self.cost = 0
        else:
            self.set_cost()

    def set_parent(self, new_parent):
        self.parent = new_parent
        self.set_cost()

    def set_cost(self):
        step_cost = sqrt((self.pos[0]-self.parent.pos[0])**2 + (self.pos[1]-self.parent.pos[1])**2)
        if world_collisions.is_seating_location(self.pos):
            step_cost *= AGENT_SEATING_PENALTY
        self.cost = self.parent.cost + step_cost

    def __eq__(self, other) -> bool:
        if type(other) == DijkstraNode:
            return self.pos == other.pos
        return False

    def __ne__(self, other) -> bool:
        if type(other) == DijkstraNode:
            return self.pos != other.pos
        return True

    def __lt__(self, other) -> bool:
        return self.cost < other.cost

    def __gt__(self, other) -> bool:
        return self.cost > other.cost

    def __le__(self, other) -> bool:
        return self.cost <= other.cost

    def __ge__(self, other) -> bool:
        return self.cost >= other.cost

    def __hash__(self):
        return hash(self.pos)

world_collisions = WorldManager(Image.open(WORLD_FILE))
world_collisions.read_paths()
agent_collisions = AgentCollisionManager(WORLD_EXTENTS[0], WORLD_EXTENTS[1])

def run_sim(percent_filled, time, draw_interval = 1, pathing = "random"):
    agents = []

    for gridpos in world_collisions.seating_pixels:
        for pos in [[gridpos[0]+0.25, gridpos[1]+0.25],[gridpos[0]+0.25, gridpos[1]+0.75],[gridpos[0]+0.75, gridpos[1]+0.25],[gridpos[0]+0.75, gridpos[1]+0.75]]:
            if random() < percent_filled:
                a = Agent((pos[0], pos[1]), AGENT_STEP, world_collisions.assign_goal(pathing,pos), 0)
                agents.append(a)
                agent_collisions.register_member(a)
    print(f"simulating {len(agents)} agents")

    agent_positions = [[] for _ in range(len(agents))]
    agent_counts = []

    plt.imshow(plt.imread(WORLD_FILE), extent=[0, WORLD_EXTENTS[0], 0, WORLD_EXTENTS[1]])
    plt.savefig("temp.png")
    plt.clf()
    frame = cv2.imread("temp.png")
    height, width, _ = frame.shape

    fourcc = cv2.VideoWriter_fourcc(*'mp4v')
    video = cv2.VideoWriter("movement.mp4", fourcc, 30, (width, height))

    for t in range(time):
        print(str((t+1)/time * 100) + "%")
        done = True
        c = 0
        for i, a in enumerate(agents):
            if not a.finished:
                c += 1
                done = False
                a.update()
                if a.finished:
                    agent_positions[i].append([0,0])
                    continue
                if t%draw_interval == 0:
                    agent_positions[i].append((a.x_pos, a.y_pos))
        agent_counts.append(c)

        if done:
            break

        if t%draw_interval == 0:
            plt.scatter(list(map(lambda x: x[-1][0], agent_positions)), list(map(lambda x: x[-1][1], agent_positions)), s=4)
            plt.imshow(plt.imread(WORLD_FILE), extent=[0, WORLD_EXTENTS[0], 0, WORLD_EXTENTS[1]])
            # plt.show()
            plt.savefig("temp.png")
            plt.clf()
            video.write(cv2.imread("temp.png"))

    cv2.destroyAllWindows()
    video.release()

    for pos_list in agent_positions:
        plt.plot(list(map(lambda x: x[0], pos_list[:-1])), list(map(lambda x: x[1], pos_list[:-1])), 'o-', linewidth=1, markersize=3)
    plt.imshow(plt.imread(WORLD_FILE), extent=[0, WORLD_EXTENTS[0], 0, WORLD_EXTENTS[1]])
    plt.show()

    print(agent_counts)

    return (t+1) * AGENT_STEP/1.4

if __name__ == "__main__":
    max_time = 30*60*4 # ~half an hour

    # cProfile.run("run_sim(1.0, max_time, draw_interval = 20)","profilestats")
    with open("results.txt","w") as f:
        random_100 = run_sim(1.0, max_time, draw_interval = 10, pathing = 'random')
        print(f"the agents exited the stadium in {random_100} seconds")
        f.write(str(random_100) + "\n")


text:
"""
TODO: Test insert_data under several modes
"""

# pylint: disable=redefined-outer-name, unused-argument

from __future__ import division, print_function
import os
import tempfile
import uuid
import string

import numpy as np
from osgeo import gdal
import pytest
import shapely.ops

import buzzard as buzz
from buzzard.test.tools import fpeq, sreq, eq
from buzzard.test import make_tile_set
from .tools import  get_srs_by_name

SR1 = get_srs_by_name('EPSG:2154')
SR2 = get_srs_by_name('EPSG:2154 centered')

# FIXTURES ************************************************************************************** **
@pytest.fixture(scope='module')
def fps():
    """See make_tile_set
    A B C
    D E F
    G H I
    """
    return make_tile_set.make_tile_set(3, [0.1, -0.1])

@pytest.fixture()
def env():
    with buzz.Env(significant=8):
        yield

@pytest.fixture(scope='module')
def shp1_path(fps):
    """Create a shapefile in SR1 containing all single letter polygons from `fps` fixture"""
    path = '{}/{}.shp'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    ds.create_vector('poly', path, 'polygon', sr=SR1['wkt'])
    for letter in string.ascii_uppercase[:9]:
        ds.poly.insert_data(fps[letter].poly)
    ds.poly.close()
    del ds
    yield path
    gdal.GetDriverByName('ESRI Shapefile').Delete(path)

@pytest.fixture(scope='module')
def tif1_path(fps):
    """Create a tif in SR1 with all single letter footprints from `fps` fixture burnt in it"""
    path = '{}/{}.tif'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    ds.create_raster('rast', path, fps.AI, 'int32', 1, sr=SR1['wkt'])
    for letter in string.ascii_uppercase[:9]:
        fp = fps[letter]
        arr = np.full(fp.shape, ord(letter), dtype=int)
        ds.rast.set_data(arr, fp=fp)
    ds.rast.close()
    del ds
    yield path
    gdal.GetDriverByName('GTiff').Delete(path)

@pytest.fixture(scope='module')
def shp2_path(fps):
    """Create a shapefile in SR2 containing all single letter polygons from `fps` fixture"""
    path = '{}/{}.shp'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    ds.create_vector('poly', path, 'polygon', sr=SR2['wkt'])
    for letter in string.ascii_uppercase[:9]:
        ds.poly.insert_data(fps[letter].poly)
    ds.poly.close()
    del ds
    yield path
    gdal.GetDriverByName('ESRI Shapefile').Delete(path)

@pytest.fixture(scope='module')
def tif2_path(fps):
    """Create a tif in SR2 with all single letter footprints from `fps` fixture burnt in it"""
    path = '{}/{}.tif'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    with ds.acreate_raster(path, fps.AI, 'int32', 1, sr=SR2['wkt']).close as r:
        for letter in string.ascii_uppercase[:9]:
            fp = fps[letter]
            arr = np.full(fp.shape, ord(letter), dtype=int)
            r.set_data(arr, fp=fp)
    yield path
    gdal.GetDriverByName('GTiff').Delete(path)

@pytest.fixture(scope='module')
def shp3_path(fps):
    """Create a shapefile without SR containing all single letter polygons from `fps` fixture"""
    path = '{}/{}.shp'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    ds.create_vector('poly', path, 'polygon', sr=None)
    for letter in string.ascii_uppercase[:9]:
        ds.poly.insert_data(fps[letter].poly)
    ds.poly.close()
    del ds
    yield path
    gdal.GetDriverByName('ESRI Shapefile').Delete(path)

@pytest.fixture(scope='module')
def tif3_path(fps):
    """Create a tif witwhout SR with all single letter footprints from `fps` fixture burnt in it"""
    path = '{}/{}.tif'.format(tempfile.gettempdir(), uuid.uuid4())

    ds = buzz.DataSource()
    with ds.acreate_raster(path, fps.AI, 'int32', 1, sr=None).close as r:
        for letter in string.ascii_uppercase[:9]:
            fp = fps[letter]
            arr = np.full(fp.shape, ord(letter), dtype=int)
            r.set_data(arr, fp=fp)
    yield path
    gdal.GetDriverByName('GTiff').Delete(path)

@pytest.fixture()
def random_path_shp():
    """Create a temporary path, and take care of cleanup afterward"""
    path = '{}/{}.shp'.format(tempfile.gettempdir(), uuid.uuid4())
    yield path
    if os.path.isfile(path):
        try:
            gdal.GetDriverByName('ESRI Shapefile').Delete(path)
        except:
            os.remove(path)

@pytest.fixture()
def random_path_tif():
    """Create a temporary path, and take care of cleanup afterward"""
    path = '{}/{}.tif'.format(tempfile.gettempdir(), uuid.uuid4())
    yield path
    if os.path.isfile(path):
        try:
            gdal.GetDriverByName('GTiff').Delete(path)
        except:
            os.remove(path)

# TESTS ***************************************************************************************** **
def test_mode1(fps, shp1_path, tif1_path, shp2_path, tif2_path, shp3_path, tif3_path):
    ds = buzz.DataSource()
    ds.open_raster('tif1', tif1_path)
    ds.open_vector('shp1', shp1_path)
    ds.open_raster('tif2', tif2_path)
    ds.open_vector('shp2', shp2_path)
    ds.open_raster('tif3', tif3_path)
    ds.open_vector('shp3', shp3_path)

    # Test SR equality
    assert sreq(
        ds.tif1.wkt_virtual,
        ds.tif1.wkt_stored,
        ds.tif1.proj4_virtual,
        ds.tif1.proj4_stored,
        ds.shp1.wkt_virtual,
        ds.shp1.wkt_stored,
        ds.shp1.proj4_virtual,
        ds.shp1.proj4_stored,
    )
    assert sreq(
        ds.tif2.wkt_virtual,
        ds.tif2.wkt_stored,
        ds.tif2.proj4_virtual,
        ds.tif2.proj4_stored,
        ds.shp2.wkt_virtual,
        ds.shp2.wkt_stored,
        ds.shp2.proj4_virtual,
        ds.shp2.proj4_stored,
    )
    assert not sreq(ds.tif1.wkt_stored, ds.tif2.wkt_stored)
    assert (
        None ==
        ds.wkt ==
        ds.proj4 ==
        ds.tif3.wkt_virtual ==
        ds.tif3.wkt_stored ==
        ds.tif3.proj4_virtual ==
        ds.tif3.proj4_stored ==
        ds.shp3.wkt_virtual ==
        ds.shp3.wkt_stored ==
        ds.shp3.proj4_virtual ==
        ds.shp3.proj4_stored
    )

    # Test footprints equality
    assert fpeq(
        fps.AI,
        # tif/shp 1
        ds.tif1.fp,
        ds.tif1.fp_origin,
        buzz.Footprint.of_extent(ds.shp1.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds[[0, 2, 1, 3]], fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
        # tif/shp 2
        ds.tif2.fp,
        ds.tif2.fp_origin,
        buzz.Footprint.of_extent(ds.shp2.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds[[0, 2, 1, 3]], fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
        # tif/shp 3
        ds.tif3.fp,
        ds.tif3.fp_origin,
        buzz.Footprint.of_extent(ds.shp3.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds[[0, 2, 1, 3]], fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
    )

    # Test what's written in all 6 files
    tif1 = ds.tif1.get_data()
    tif2 = ds.tif2.get_data()
    tif3 = ds.tif3.get_data()
    assert np.all(tif1 == tif2)
    assert np.all(tif1 == tif3)
    for i, letter in enumerate(string.ascii_uppercase[:9]):
        # tif/shp 1
        shp1 = ds.shp1.get_data(i, None)
        raster1_polys = ds.tif1.fp.find_polygons(tif1 == ord(letter))
        assert len(raster1_polys) == 1
        assert (shp1 ^ raster1_polys[0]).is_empty

        # tif/shp 2
        shp2 = ds.shp2.get_data(i, None)
        raster2_polys = ds.tif2.fp.find_polygons(tif2 == ord(letter))
        assert len(raster2_polys) == 1
        assert (shp2 ^ raster2_polys[0]).is_empty

        # tif/shp 3
        shp3 = ds.shp3.get_data(i, None)
        raster3_polys = ds.tif3.fp.find_polygons(tif3 == ord(letter))
        assert len(raster3_polys) == 1
        assert (shp3 ^ raster3_polys[0]).is_empty

def test_mode2(fps, shp1_path, tif1_path, shp2_path, tif2_path, shp3_path, tif3_path,
               random_path_shp, random_path_tif, env):
    ds = buzz.DataSource(sr_work=SR1['wkt'])
    ds.open_raster('tif1', tif1_path)
    ds.open_vector('shp1', shp1_path)
    ds.open_raster('tif2', tif2_path)
    ds.open_vector('shp2', shp2_path)

    # Test file creation/opening without spatial reference
    with buzz.Env(allow_complex_footprint=True):
        with pytest.raises(ValueError, match='spatial refe'):
            ds.acreate_vector(random_path_shp, 'polygon', [], sr=None)
        with pytest.raises(ValueError, match='spatial refe'):
            ds.acreate_raster(random_path_tif, fps.AI, 'int32', 1, {}, sr=None)
        with pytest.raises(ValueError, match='spatial refe'):
            ds.aopen_raster(tif3_path)
        with pytest.raises(ValueError, match='spatial refe'):
            ds.aopen_vector(shp3_path)

    # Test SR equality
    assert sreq(
        ds.wkt,
        ds.proj4,
        ds.tif1.wkt_virtual,
        ds.tif1.wkt_stored,
        ds.tif1.proj4_virtual,
        ds.tif1.proj4_stored,
        ds.shp1.wkt_virtual,
        ds.shp1.wkt_stored,
        ds.shp1.proj4_virtual,
        ds.shp1.proj4_stored,
    )
    assert sreq(
        ds.tif2.wkt_virtual,
        ds.tif2.proj4_virtual,
        ds.shp2.wkt_virtual,
        ds.shp2.proj4_virtual,
        ds.tif2.wkt_stored,
        ds.tif2.proj4_stored,
        ds.shp2.wkt_stored,
        ds.shp2.proj4_stored,
    )
    assert not sreq(ds.tif1.wkt_stored, ds.tif2.wkt_stored)

    # Test foorprints equality
    assert fpeq(
        fps.AI,
        # tif/shp 1
        ds.tif1.fp,
        ds.tif1.fp_origin,
        buzz.Footprint.of_extent(ds.shp1.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds[[0, 2, 1, 3]], fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 2
        ds.tif2.fp_origin,
        buzz.Footprint.of_extent(ds.shp2.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert fpeq(
        ds.tif2.fp,
        buzz.Footprint.of_extent(ds.shp2.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert ds.tif2.fp != ds.tif1.fp

    # Test file creation with/without conversion of footprint
    with buzz.Env(allow_complex_footprint=True):
        with ds.acreate_raster(random_path_tif, fps.AI, 'int32', 1, sr=SR1['wkt']).delete as r:
            assert fpeq(
                fps.AI,
                r.fp,
                r.fp_origin
            )
        with ds.acreate_raster(random_path_tif, fps.AI, 'int32', 1, sr=SR2['wkt']).delete as r:
            assert fpeq(
                fps.AI,
                r.fp,
            )
            assert fps.AI != r.fp_origin

    # Test what's written in all 4 files
    tif1 = ds.tif1.get_data()
    tif2 = ds.tif2.get_data()
    assert np.all(tif1 == tif2)

    def f(x, y, z=None):
        return np.around(x, 6), np.around(y, 6)

    for i, letter in enumerate(string.ascii_uppercase[:9]):
        # tif/shp 1
        shp1 = ds.shp1.get_data(i, None)
        shp1 = shapely.ops.transform(f, shp1)

        raster1_polys = ds.tif1.fp.find_polygons(tif1 == ord(letter))
        assert len(raster1_polys) == 1
        raster1_poly = raster1_polys[0]
        raster1_poly = shapely.ops.transform(f, shp1)

        assert (shp1 ^ raster1_poly).is_empty

        # tif/shp 2
        shp2 = ds.shp2.get_data(i, None)
        shp2 = shapely.ops.transform(f, shp2)

        raster2_polys = ds.tif2.fp.find_polygons(tif2 == ord(letter))
        assert len(raster2_polys) == 1
        raster2_poly = raster2_polys[0]
        raster2_poly = shapely.ops.transform(f, shp2)

        assert (shp2 ^ raster2_poly).is_empty

def test_mode3(fps, shp1_path, tif1_path, shp2_path, tif2_path, shp3_path, tif3_path,
               random_path_shp, random_path_tif, env):
    ds = buzz.DataSource(sr_work=SR1['wkt'], sr_fallback=SR2['wkt'])
    ds.open_raster('tif1', tif1_path)
    ds.open_vector('shp1', shp1_path)
    ds.open_raster('tif2', tif2_path)
    ds.open_vector('shp2', shp2_path)
    ds.open_raster('tif3', tif3_path)
    ds.open_vector('shp3', shp3_path)

    # Test file creation without spatial reference
    with buzz.Env(allow_complex_footprint=True):
        with ds.acreate_vector(random_path_shp, 'polygon', [], sr=None).close as r:
            assert r.wkt_stored == None
            assert sreq(r.wkt_virtual, SR2['wkt'])
        with ds.acreate_raster(random_path_tif, fps.AI, 'int32', 1, {}, sr=None).close as v:
            assert v.wkt_stored == None
            assert sreq(v.wkt_virtual, SR2['wkt'])

    # Test SR equality
    assert sreq(
        ds.wkt,
        ds.proj4,
        ds.tif1.wkt_virtual,
        ds.tif1.wkt_stored,
        ds.tif1.proj4_virtual,
        ds.tif1.proj4_stored,
        ds.shp1.wkt_virtual,
        ds.shp1.wkt_stored,
        ds.shp1.proj4_virtual,
        ds.shp1.proj4_stored,

    )
    assert sreq(
        ds.tif2.wkt_virtual,
        ds.tif2.proj4_virtual,
        ds.shp2.wkt_virtual,
        ds.shp2.proj4_virtual,
        ds.tif2.wkt_stored,
        ds.tif2.proj4_stored,
        ds.shp2.wkt_stored,
        ds.shp2.proj4_stored,

        ds.tif3.wkt_virtual,
        ds.tif3.proj4_virtual,
        ds.shp3.wkt_virtual,
        ds.shp3.proj4_virtual,
    )
    assert not sreq(ds.tif1.wkt_stored, ds.tif2.wkt_stored)
    assert (
        None ==
        ds.tif3.wkt_stored ==
        ds.tif3.proj4_stored ==
        ds.shp3.wkt_stored ==
        ds.shp3.proj4_stored
    )

    # Test foorprints equality
    assert fpeq(
        fps.AI,
        # tif/shp 1
        ds.tif1.fp,
        ds.tif1.fp_origin,
        buzz.Footprint.of_extent(ds.shp1.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds[[0, 2, 1, 3]], fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 2
        ds.tif2.fp_origin,
        buzz.Footprint.of_extent(ds.shp2.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 3
        ds.tif3.fp_origin,
        buzz.Footprint.of_extent(ds.shp3.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert fpeq(
        ds.tif2.fp,
        buzz.Footprint.of_extent(ds.shp2.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds[[0, 2, 1, 3]], fps.AI.scale),

        ds.tif3.fp,
        buzz.Footprint.of_extent(ds.shp3.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert ds.tif2.fp != ds.tif1.fp

    # Test what's written in all 6 files
    tif1 = ds.tif1.get_data()
    tif2 = ds.tif2.get_data()
    tif3 = ds.tif3.get_data()
    assert np.all(tif1 == tif2)
    assert np.all(tif1 == tif3)

    def f(x, y, z=None):
        return np.around(x, 6), np.around(y, 6)

    for i, letter in enumerate(string.ascii_uppercase[:9]):
        # tif/shp 1
        shp1 = ds.shp1.get_data(i, None)
        shp1 = shapely.ops.transform(f, shp1)

        raster1_polys = ds.tif1.fp.find_polygons(tif1 == ord(letter))
        assert len(raster1_polys) == 1
        raster1_poly = raster1_polys[0]
        raster1_poly = shapely.ops.transform(f, shp1)

        assert (shp1 ^ raster1_poly).is_empty

        # tif/shp 2
        shp2 = ds.shp2.get_data(i, None)
        shp2 = shapely.ops.transform(f, shp2)

        raster2_polys = ds.tif2.fp.find_polygons(tif2 == ord(letter))
        assert len(raster2_polys) == 1
        raster2_poly = raster2_polys[0]
        raster2_poly = shapely.ops.transform(f, shp2)

        assert (shp2 ^ raster2_poly).is_empty

        # tif/shp 3
        shp3 = ds.shp3.get_data(i, None)
        shp3 = shapely.ops.transform(f, shp3)

        raster3_polys = ds.tif3.fp.find_polygons(tif3 == ord(letter))
        assert len(raster3_polys) == 1
        raster3_poly = raster3_polys[0]
        raster3_poly = shapely.ops.transform(f, shp3)

        assert (shp3 ^ raster3_poly).is_empty

def test_mode4(fps, shp1_path, tif1_path, shp2_path, tif2_path, shp3_path, tif3_path,
               random_path_shp, random_path_tif, env):
    ds = buzz.DataSource(sr_work=SR1['wkt'], sr_forced=SR2['wkt'])
    ds.open_raster('tif1', tif1_path)
    ds.open_vector('shp1', shp1_path)
    ds.open_raster('tif2', tif2_path)
    ds.open_vector('shp2', shp2_path)
    ds.open_raster('tif3', tif3_path)
    ds.open_vector('shp3', shp3_path)

    # Test file creation without spatial reference
    with buzz.Env(allow_complex_footprint=True):
        with ds.acreate_vector(random_path_shp, 'polygon', [], sr=None).close as r:
            assert r.wkt_stored == None
            assert sreq(r.wkt_virtual, SR2['wkt'])
        with ds.acreate_raster(random_path_tif, fps.AI, 'int32', 1, {}, sr=None).close as v:
            assert v.wkt_stored == None
            assert sreq(v.wkt_virtual, SR2['wkt'])

    # Test SR equality
    assert sreq(
        ds.wkt,
        ds.proj4,
        ds.tif1.wkt_stored,
        ds.tif1.proj4_stored,
        ds.shp1.wkt_stored,
        ds.shp1.proj4_stored,
    )
    assert sreq(
        ds.tif1.wkt_virtual,
        ds.tif1.proj4_virtual,
        ds.shp1.wkt_virtual,
        ds.shp1.proj4_virtual,

        ds.tif2.wkt_virtual,
        ds.tif2.proj4_virtual,
        ds.shp2.wkt_virtual,
        ds.shp2.proj4_virtual,
        ds.tif2.wkt_stored,
        ds.tif2.proj4_stored,
        ds.shp2.wkt_stored,
        ds.shp2.proj4_stored,

        ds.tif3.wkt_virtual,
        ds.tif3.proj4_virtual,
        ds.shp3.wkt_virtual,
        ds.shp3.proj4_virtual,
    )
    assert not sreq(ds.tif1.wkt_stored, ds.tif2.wkt_stored)
    assert (
        None ==
        ds.tif3.wkt_stored ==
        ds.tif3.proj4_stored ==
        ds.shp3.wkt_stored ==
        ds.shp3.proj4_stored
    )

    # Test foorprints equality
    assert fpeq(
        fps.AI,
        # tif/shp 1
        ds.tif1.fp_origin,
        buzz.Footprint.of_extent(ds.shp1.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 2
        ds.tif2.fp_origin,
        buzz.Footprint.of_extent(ds.shp2.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 3
        ds.tif3.fp_origin,
        buzz.Footprint.of_extent(ds.shp3.extent_stored, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds_stored[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert fpeq(
        # tif/shp 1
        ds.tif1.fp,
        buzz.Footprint.of_extent(ds.shp1.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp1.bounds[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 2
        ds.tif2.fp,
        buzz.Footprint.of_extent(ds.shp2.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp2.bounds[[0, 2, 1, 3]], fps.AI.scale),

        # tif/shp 3
        ds.tif3.fp,
        buzz.Footprint.of_extent(ds.shp3.extent, fps.AI.scale),
        buzz.Footprint.of_extent(ds.shp3.bounds[[0, 2, 1, 3]], fps.AI.scale),
    )
    assert ds.tif1.fp != ds.tif1.fp_stored

    # Test what's written in all 6 files
    tif1 = ds.tif1.get_data()
    tif2 = ds.tif2.get_data()
    tif3 = ds.tif3.get_data()
    assert np.all(tif1 == tif2)
    assert np.all(tif1 == tif3)

    def f(x, y, z=None):
        return np.around(x, 6), np.around(y, 6)

    for i, letter in enumerate(string.ascii_uppercase[:9]):
        # tif/shp 1
        shp1 = ds.shp1.get_data(i, None)
        shp1 = shapely.ops.transform(f, shp1)

        raster1_polys = ds.tif1.fp.find_polygons(tif1 == ord(letter))
        assert len(raster1_polys) == 1
        raster1_poly = raster1_polys[0]
        raster1_poly = shapely.ops.transform(f, shp1)

        assert (shp1 ^ raster1_poly).is_empty

        # tif/shp 2
        shp2 = ds.shp2.get_data(i, None)
        shp2 = shapely.ops.transform(f, shp2)

        raster2_polys = ds.tif2.fp.find_polygons(tif2 == ord(letter))
        assert len(raster2_polys) == 1
        raster2_poly = raster2_polys[0]
        raster2_poly = shapely.ops.transform(f, shp2)

        assert (shp2 ^ raster2_poly).is_empty

        # tif/shp 3
        shp3 = ds.shp3.get_data(i, None)
        shp3 = shapely.ops.transform(f, shp3)

        raster3_polys = ds.tif3.fp.find_polygons(tif3 == ord(letter))
        assert len(raster3_polys) == 1
        raster3_poly = raster3_polys[0]
        raster3_poly = shapely.ops.transform(f, shp3)

        assert (shp3 ^ raster3_poly).is_empty

max_stars_repo_name: azarz/buzzard

max_issues_repo_name: azarz/buzzard

max_forks_repo_name: azarz/buzzard


text:
"""Longley dataset"""

__all__ = ['COPYRIGHT','TITLE','SOURCE','DESCRSHORT','DESCRLONG','NOTE', 'load']

__docformat__ = 'restructuredtext'

COPYRIGHT   = """This is public domain."""
TITLE       = __doc__
SOURCE      = """
The classic 1967 Longley Data

http://www.itl.nist.gov/div898/strd/lls/data/Longley.shtml

::

    Longley, J.W. (1967) "An Appraisal of Least Squares Programs for the
        Electronic Comptuer from the Point of View of the User."  Journal of
        the American Statistical Association.  62.319, 819-41.
"""

DESCRSHORT  = """"""

DESCRLONG   = """The Longley dataset contains various US macroeconomic
variables that are known to be highly collinear.  It has been used to appraise
the accuracy of least squares routines."""

NOTE        = """
Number of Observations - 16

Number of Variables - 6

Variable name definitions::

        TOTEMP - Total Employment
        GNPDEFL - GNP deflator
        GNP - GNP
        UNEMP - Number of unemployed
        ARMED - Size of armed forces
        POP - Population
        YEAR - Year (1947 - 1962)
"""

from numpy import recfromtxt, array, column_stack
from scikits.statsmodels.datasets import Dataset
from os.path import dirname, abspath

def load():
    """
    Load the Longley data and return a Dataset class.

    Returns
    -------
    Dataset instance
        See DATASET_PROPOSAL.txt for more information.
    """
    filepath = dirname(abspath(__file__))
    data = recfromtxt(open(filepath+'/longley.csv',"rb"), delimiter=",", names=True,
            dtype=float, usecols=(1,2,3,4,5,6,7))
    names = list(data.dtype.names)
    endog = array(data[names[0]], dtype=float)
    endog_name = names[0]
    exog = column_stack(data[i] for i in names[1:]).astype(float)
    exog_name = names[1:]
    dataset = Dataset(data=data, names=names, endog=endog, exog=exog,
            endog_name=endog_name, exog_name = exog_name)
    return dataset

max_stars_repo_name: matthew-brett/statsmodels

max_issues_repo_name: matthew-brett/statsmodels

max_forks_repo_name: matthew-brett/statsmodels

reason: from numpy