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

# This code is part of Qiskit.
#
# (C) Copyright IBM 2018, 2020.
#
# This code is licensed under the Apache License, Version 2.0. You may
# obtain a copy of this license in the LICENSE.txt file in the root directory
# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
#
# Any modifications or derivative works of this code must retain this
# copyright notice, and modified files need to carry a notice indicating
# that they have been altered from the originals.
"""
The Amplitude Estimation Algorithm.
"""

from typing import Optional
import logging
from collections import OrderedDict
import numpy as np
from scipy.stats import chi2, norm
from scipy.optimize import bisect

from qiskit.aqua import AquaError
from qiskit.aqua.utils import CircuitFactory
from qiskit.aqua.circuits import PhaseEstimationCircuit
from qiskit.aqua.components.iqfts import IQFT, Standard
from qiskit.aqua.utils.validation import validate_min
from .ae_algorithm import AmplitudeEstimationAlgorithm
from .ae_utils import pdf_a, derivative_log_pdf_a, bisect_max

logger = logging.getLogger(__name__)

# pylint: disable=invalid-name

class AmplitudeEstimation(AmplitudeEstimationAlgorithm):
    """
    The Amplitude Estimation algorithm.
    """

    def __init__(self, num_eval_qubits: int,
                 a_factory: Optional[CircuitFactory] = None,
                 i_objective: Optional[int] = None,
                 q_factory: Optional[CircuitFactory] = None,
                 iqft: Optional[IQFT] = None) -> None:
        """

        Args:
            num_eval_qubits: number of evaluation qubits, has a min. value of 1.
            a_factory: the CircuitFactory subclass object representing
                                        the problem unitary
            i_objective: i objective
            q_factory: the CircuitFactory subclass object representing an
                                        amplitude estimation sample (based on a_factory)
            iqft: the Inverse Quantum Fourier Transform component,
                         defaults to using a standard iqft when None
        """
        validate_min('num_eval_qubits', num_eval_qubits, 1)
        super().__init__(a_factory, q_factory, i_objective)

        # get parameters
        self._m = num_eval_qubits
        self._M = 2 ** num_eval_qubits

        if iqft is None:
            iqft = Standard(self._m)

        self._iqft = iqft
        self._circuit = None
        self._ret = {}

    @property
    def _num_qubits(self):
        if self.a_factory is None:  # if A factory is not set, no qubits are specified
            return 0

        num_ancillas = self.q_factory.required_ancillas_controlled()
        num_qubits = self.a_factory.num_target_qubits + self._m + num_ancillas

        return num_qubits

    def construct_circuit(self, measurement=False):
        """
        Construct the Amplitude Estimation quantum circuit.

        Args:
            measurement (bool): Boolean flag to indicate if measurement
                should be included in the circuit.

        Returns:
            QuantumCircuit: the QuantumCircuit object for the constructed circuit
        """
        pec = PhaseEstimationCircuit(
            iqft=self._iqft, num_ancillae=self._m,
            state_in_circuit_factory=self.a_factory,
            unitary_circuit_factory=self.q_factory
        )

        self._circuit = pec.construct_circuit(measurement=measurement)
        return self._circuit

    def _evaluate_statevector_results(self, probabilities):
        # map measured results to estimates
        y_probabilities = OrderedDict()
        for i, probability in enumerate(probabilities):
            b = '{0:b}'.format(i).rjust(self._num_qubits, '0')[::-1]
            y = int(b[:self._m], 2)
            y_probabilities[y] = y_probabilities.get(y, 0) + probability

        a_probabilities = OrderedDict()
        for y, probability in y_probabilities.items():
            if y >= int(self._M / 2):
                y = self._M - y
            a = np.round(np.power(np.sin(y * np.pi / 2 ** self._m), 2),
                         decimals=7)
            a_probabilities[a] = a_probabilities.get(a, 0) + probability

        return a_probabilities, y_probabilities

    def _compute_fisher_information(self, observed=False):
        fisher_information = None
        mlv = self._ret['ml_value']  # MLE in [0,1]
        m = self._m
        if observed:
            ai = np.asarray(self._ret['values'])
            pi = np.asarray(self._ret['probabilities'])

            # Calculate the observed Fisher information
            fisher_information = sum(p * derivative_log_pdf_a(a, mlv, m)**2 for p, a in zip(pi, ai))
        else:
            def integrand(x):
                return (derivative_log_pdf_a(x, mlv, m))**2 * pdf_a(x, mlv, m)

            M = 2**m
            grid = np.sin(np.pi * np.arange(M / 2 + 1) / M)**2
            fisher_information = sum(integrand(x) for x in grid)

        return fisher_information

    def _fisher_ci(self, alpha, observed=False):
        shots = self._ret['shots']
        mle = self._ret['ml_value']

        std = np.sqrt(shots * self._compute_fisher_information(observed))
        ci = mle + norm.ppf(1 - alpha / 2) / std * np.array([-1, 1])

        return [self.a_factory.value_to_estimation(bound) for bound in ci]

    def _likelihood_ratio_ci(self, alpha):
        # Compute the two intervals in which we the look for values above
        # the likelihood ratio: the two bubbles next to the QAE estimate
        M = 2**self._m
        qae = self._ret['value']
        y = M * np.arcsin(np.sqrt(qae)) / np.pi
        left_of_qae = np.sin(np.pi * (y - 1) / M)**2
        right_of_qae = np.sin(np.pi * (y + 1) / M)**2

        bubbles = [left_of_qae, qae, right_of_qae]

        # likelihood function
        ai = np.asarray(self._ret['values'])
        pi = np.asarray(self._ret['probabilities'])
        m = self._m
        shots = self._ret['shots']

        def loglikelihood(a):
            return np.sum(shots * pi * np.log(pdf_a(ai, a, m)))

        # The threshold above which the likelihoods are in the
        # confidence interval
        loglik_mle = loglikelihood(self._ret['ml_value'])
        thres = loglik_mle - chi2.ppf(1 - alpha, df=1) / 2

        def cut(x):
            return loglikelihood(x) - thres

        # Store the boundaries of the confidence interval
        lower = upper = self._ret['ml_value']

        # Check the two intervals/bubbles: check if they surpass the
        # threshold and if yes add the part that does to the CI
        for a, b in zip(bubbles[:-1], bubbles[1:]):
            # Compute local maximum and perform a bisect search between
            # the local maximum and the bubble boundaries
            locmax, val = bisect_max(loglikelihood, a, b, retval=True)
            if val >= thres:
                # Bisect pre-condition is that the function has different
                # signs at the boundaries of the interval we search in
                if cut(a) * cut(locmax) < 0:
                    left = bisect(cut, a, locmax)
                    lower = np.minimum(lower, left)
                if cut(locmax) * cut(b) < 0:
                    right = bisect(cut, locmax, b)
                    upper = np.maximum(upper, right)

        # Put together CI
        ci = [lower, upper]
        return [self.a_factory.value_to_estimation(bound) for bound in ci]

    def confidence_interval(self, alpha, kind='likelihood_ratio'):
        """
        Compute the (1 - alpha) confidence interval

        Args:
            alpha (float): confidence level: compute the (1 - alpha) confidence interval
            kind (str): the method to compute the confidence interval, can be 'fisher',
                'observed_fisher' or 'likelihood_ratio' (default)

        Returns:
            list[float]: the (1 - alpha) confidence interval

        Raises:
            AquaError: if 'mle' is not in self._ret.keys() (i.e. `run` was not called yet)
            NotImplementedError: if the confidence interval method `kind` is not implemented
        """
        # check if AE did run already
        if 'mle' not in self._ret.keys():
            raise AquaError('Call run() first!')

        # if statevector simulator the estimate is exact
        if self._quantum_instance.is_statevector:
            return 2 * [self._ret['estimation']]

        if kind in ['likelihood_ratio', 'lr']:
            return self._likelihood_ratio_ci(alpha)

        if kind in ['fisher', 'fi']:
            return self._fisher_ci(alpha, observed=False)

        if kind in ['observed_fisher', 'observed_information', 'oi']:
            return self._fisher_ci(alpha, observed=True)

        raise NotImplementedError('CI `{}` is not implemented.'.format(kind))

    def _run_mle(self):
        """
        Compute the Maximum Likelihood Estimator (MLE)

        Returns:
            The MLE for the previous AE run

        Note: Before calling this method, call the method `run` of the AmplitudeEstimation instance
        """
        M = self._M
        qae = self._ret['value']

        # likelihood function
        ai = np.asarray(self._ret['values'])
        pi = np.asarray(self._ret['probabilities'])
        m = self._m
        shots = self._ret['shots']

        def loglikelihood(a):
            return np.sum(shots * pi * np.log(pdf_a(ai, a, m)))

        # y is pretty much an integer, but to map 1.9999 to 2 we must first
        # use round and then int conversion
        y = int(np.round(M * np.arcsin(np.sqrt(qae)) / np.pi))

        # Compute the two intervals in which are candidates for containing
        # the maximum of the log-likelihood function: the two bubbles next to
        # the QAE estimate
        bubbles = None
        if y == 0:
            right_of_qae = np.sin(np.pi * (y + 1) / M)**2
            bubbles = [qae, right_of_qae]

        elif y == int(M / 2):
            left_of_qae = np.sin(np.pi * (y - 1) / M)**2
            bubbles = [left_of_qae, qae]

        else:
            left_of_qae = np.sin(np.pi * (y - 1) / M)**2
            right_of_qae = np.sin(np.pi * (y + 1) / M)**2
            bubbles = [left_of_qae, qae, right_of_qae]

        # Find global maximum amongst the two local maxima
        a_opt = qae
        loglik_opt = loglikelihood(a_opt)
        for a, b in zip(bubbles[:-1], bubbles[1:]):
            locmax, val = bisect_max(loglikelihood, a, b, retval=True)
            if val > loglik_opt:
                a_opt = locmax
                loglik_opt = val

        # Convert the value to an estimation
        val_opt = self.a_factory.value_to_estimation(a_opt)

        # Store MLE and the MLE mapped to an estimation
        self._ret['ml_value'] = a_opt
        self._ret['mle'] = val_opt

    def _run(self):
        # check if A factory has been set
        if self.a_factory is None:
            raise AquaError("a_factory must be set!")

        if self._quantum_instance.is_statevector:
            self.construct_circuit(measurement=False)
            # run circuit on statevector simulator
            ret = self._quantum_instance.execute(self._circuit)
            state_vector = np.asarray([ret.get_statevector(self._circuit)])
            self._ret['statevector'] = state_vector

            # get state probabilities
            state_probabilities = np.real(
                state_vector.conj() * state_vector)[0]

            # evaluate results
            a_probabilities, y_probabilities = self._evaluate_statevector_results(
                state_probabilities)

            # store number of shots: convention is 1 shot for statevector,
            # needed so that MLE works!
            self._ret['shots'] = 1
        else:
            # run circuit on QASM simulator
            self.construct_circuit(measurement=True)
            ret = self._quantum_instance.execute(self._circuit)

            # get counts
            self._ret['counts'] = ret.get_counts()

            # construct probabilities
            y_probabilities = {}
            a_probabilities = {}
            shots = sum(ret.get_counts().values())
            for state, counts in ret.get_counts().items():
                y = int(state.replace(' ', '')[:self._m][::-1], 2)
                p = counts / shots
                y_probabilities[y] = p
                a = np.round(np.power(np.sin(y * np.pi / 2 ** self._m), 2),
                             decimals=7)
                a_probabilities[a] = a_probabilities.get(a, 0.0) + p

            # store shots
            self._ret['shots'] = shots

        # construct a_items and y_items
        a_items = [(a, p) for (a, p) in a_probabilities.items() if p > 1e-6]
        y_items = [(y, p) for (y, p) in y_probabilities.items() if p > 1e-6]
        a_items = list(a_probabilities.items())
        y_items = list(y_probabilities.items())
        a_items = sorted(a_items)
        y_items = sorted(y_items)
        self._ret['a_items'] = a_items
        self._ret['y_items'] = y_items

        # map estimated values to original range and extract probabilities
        self._ret['mapped_values'] = [self.a_factory.value_to_estimation(
            a_item[0]) for a_item in self._ret['a_items']]
        self._ret['values'] = [a_item[0] for a_item in self._ret['a_items']]
        self._ret['y_values'] = [y_item[0] for y_item in y_items]
        self._ret['probabilities'] = [a_item[1]
                                      for a_item in self._ret['a_items']]
        self._ret['mapped_items'] = [(self._ret['mapped_values'][i], self._ret['probabilities'][i])
                                     for i in range(len(self._ret['mapped_values']))]

        # determine most likely estimator
        self._ret['value'] = None  # estimate in [0,1]
        self._ret['estimation'] = None  # estimate mapped to right interval
        self._ret['max_probability'] = 0
        for val, (est, prob) in zip(self._ret['values'], self._ret['mapped_items']):
            if prob > self._ret['max_probability']:
                self._ret['max_probability'] = prob
                self._ret['estimation'] = est
                self._ret['value'] = val

        # get MLE
        self._run_mle()

        # get 95% confidence interval
        alpha = 0.05
        kind = 'likelihood_ratio'  # empirically the most precise kind
        self._ret['95%_confidence_interval'] = self.confidence_interval(alpha, kind)

        return self._ret

max_stars_repo_name: gillenhaalb/qiskit-aqua

max_issues_repo_name: gillenhaalb/qiskit-aqua

max_forks_repo_name: gillenhaalb/qiskit-aqua

reason: import numpy,from scipy


text:
# -*- coding: utf-8 -*-
"""
Created on Mon Aug 31 04:43:44 2020

@author: Srikar Balusu
"""

import json
import pandas as pd
import numpy as np
import re
import random
import sqlite3
import datetime
import calendar
from dateutil.relativedelta import *

with open('lookup1.json') as json_file:
    data = json.load(json_file)

with open('uniquelookup.json') as json_file:
    data2 = json.load(json_file)

with open('state_dict.json') as json_file:
    state_dict = json.load(json_file)
conn = sqlite3.connect('testQ.db')
c = conn.cursor()
question_template = "What percentage of (Case Entity) in (State Entity) are from (Race Entity)?"
question_template_id = 'db2q6'
output = {}
entities = ['Case Entity', 'State Entity', 'Race Entity']
case_entity_list = ['confirmed cases', 'cases', 'deaths']
location_entity_list = ['(State)', "(State Abbreviation)", 'the state of (State)']
count = 1
question_key = {}

def genRaceInput(specific_race):
    if specific_race == 'African-Americans' or specific_race == 'Black':
        race_input = 'Black'
    elif specific_race == 'Asian':
        race_input = 'Asian'
    elif specific_race == 'White' or specific_race == 'Caucasian':
        race_input = 'White'
    elif specific_race == 'American Indian' or specific_race == 'Alaska Native' or specific_race == 'American Indian or Alaska Native':
        race_input = 'AIAN'
    elif specific_race == 'Pacific Islander' or specific_race == 'Native Hawaiian' or specific_race == 'Pacific Islander and Native Hawaiian':
        race_input = 'NHPI'
    elif specific_race == 'multiracial' or specific_race == 'mixed':
        race_input = 'Multiracial'
    else:
        race_input = 'LatinX'
    return race_input
def raceQuery(query, race_entity):
    output = race_entity
    if race_entity == 'people of color':
        race_input = 'Black'
    elif race_entity == 'people of mixed color' or race_entity == 'people of two or more races':
        race_input = 'Multiracial'
    elif race_entity.find('Non-Hispanic') >=0:
        specific_race = random.choice(data2['Race'])
        while specific_race == 'multiracial' or specific_race == 'mixed' or specific_race == 'Latino' or specific_race == 'Hispanic':
            specific_race = random.choice(data2['Race'])
        race_input = genRaceInput(specific_race)
        if specific_race == 'Pacific Islander and Native Hawaiian':
            output = output.replace("(race)", 'Pacific Islanders and Native Hawaiians')
        elif specific_race == 'American Indian or Alaska Native':
            output = output.replace("(race)", 'American Indians or Alaska Natives')
        else:
            output = output.replace("(race)", specific_race + 's')
    elif race_entity == '(race) only':
        specific_race = random.choice(data2['Race'])
        race_input = genRaceInput(specific_race)
        if specific_race == 'Pacific Islander and Native Hawaiian':
            output = output.replace("(race)", 'Pacific Islanders and Native Hawaiians')
        elif specific_race == 'American Indian or Alaska Native':
            output = output.replace("(race)", 'American Indians or Alaska Natives')
        else:
            output = output.replace("(race)", specific_race + 's')
    elif race_entity == '(race) people':
        specific_race = random.choice(data2['Race'])
        race_input = genRaceInput(specific_race)
        output = output.replace("(race)", specific_race)
    elif race_entity == '(race)':
        specific_race = random.choice(data2['Race'])
        race_input = genRaceInput(specific_race)
        if specific_race == 'Pacific Islander and Native Hawaiian':
            output = output.replace("(race)", 'Pacific Islanders and Native Hawaiians')
        elif specific_race == 'American Indian or Alaska Native':
            output = output.replace("(race)", 'American Indians or Alaska Natives')
        else:
            output = output.replace("(race)", specific_race + 's')
    query = query.replace("Race Entity Column",race_input)
    return query, output
while count < 300:
    output[count] = []
    populated_entities = []
    location = random.choice(location_entity_list)
    sql_template = "Select Case Entity Column_Race Entity Column from db2race where date = 'given date' and state = \"State Entity\" "
    race_entity = random.choice(data['Race'])
    case_entity = random.choice(case_entity_list)
    while race_entity == 'mixed' or race_entity =='multiracial':
        race_entity  = random.choice(data['Race'])

    if location.find("(State)") >= 0:
        state_name = random.choice(data2['State'])
        while state_name == 'Diamond Princess':
            state_name = random.choice(data2['State'])
        key_list = list(state_dict.keys())
        val_list = list(state_dict.values())
        state_abbreviation = key_list[val_list.index(state_name)]
        query = sql_template.replace("State Entity", state_abbreviation)
        query, race_e = raceQuery(query, race_entity)
        loc = location.replace("(State)", state_name)
    else:
        state_abbreviation = random.choice(data2['State Abbreviation'])
        query = sql_template.replace("State Entity", state_abbreviation)
        query, race_e = raceQuery(query,race_entity)
        loc = location.replace("(State Abbreviation)", state_abbreviation)
    if case_entity.find("cases") >=0:
        query = query.replace("Case Entity Column", 'Cases')
    else:
        query = query.replace("Case Entity Column", "Deaths")
    today = datetime.date.today()-datetime.timedelta(days=2)
    if today.weekday() == 1 or today.weekday() == 0:
        date = today-datetime.timedelta(days=today.weekday()+1)
    elif today.weekday() == 6 or today.weekday() == 2:
        date = today
    else:
        date = today - datetime.timedelta(days=today.weekday()-2)
    query = query.replace('given date', str(date).replace("-", ""))
    real_question = question_template.replace("(State Entity)", loc)
    real_question = real_question.replace("(Case Entity)", case_entity)
    real_question = real_question.replace("(Race Entity)", race_e)
    populated_entities.append(case_entity)
    populated_entities.append(loc)
    populated_entities.append(race_e)
    c.execute(query)
    result = c.fetchall()
    #if len(result) == 0 or result[0][0] == None:
    #    continue
    if real_question in question_key.keys():
        continue
    else:
        question_key[real_question] = True
        output[count].append({'question_template_id' : question_template_id, 'question_template' : question_template,
                     'entities' : entities, 'question' : real_question,
                 'populated_entities': populated_entities, 'query_template' : sql_template, 'query' :  query, 'database': 'database 2'})
        print(count)
        print(question_template)
        print(sql_template)
        print(real_question)
        print(query)
        print(result)

        count = count +1
with open('db2q6data.json', 'w') as outfile:
    json.dump(output,outfile)
print("done")

ext: py

max_stars_repo_name: wangpinggl/covidQA

max_issues_repo_name: wangpinggl/covidQA

max_forks_repo_name: wangpinggl/covidQA

reason: import numpy


text:
# emacs: -*- mode: python-mode; py-indent-offset: 4; indent-tabs-mode: nil -*-
# vi: set ft=python sts=4 ts=4 sw=4 et:
### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ##
#
#   See COPYING file distributed along with the NiBabel package for the
#   copyright and license terms.
#
### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ### ##
""" Processor functions for images """
import numpy as np

from .orientations import io_orientation, OrientationError
from .loadsave import load

def squeeze_image(img):
    """ Return image, remove axes length 1 at end of image shape

    For example, an image may have shape (10,20,30,1,1).  In this case
    squeeze will result in an image with shape (10,20,30).  See doctests
    for further description of behavior.

    Parameters
    ----------
    img : ``SpatialImage``

    Returns
    -------
    squeezed_img : ``SpatialImage``
       Copy of img, such that data, and data shape have been squeezed,
       for dimensions > 3rd, and at the end of the shape list

    Examples
    --------
    >>> import nibabel as nf
    >>> shape = (10,20,30,1,1)
    >>> data = np.arange(np.prod(shape)).reshape(shape)
    >>> affine = np.eye(4)
    >>> img = nf.Nifti1Image(data, affine)
    >>> img.shape == (10, 20, 30, 1, 1)
    True
    >>> img2 = squeeze_image(img)
    >>> img2.shape == (10, 20, 30)
    True

    If the data are 3D then last dimensions of 1 are ignored

    >>> shape = (10,1,1)
    >>> data = np.arange(np.prod(shape)).reshape(shape)
    >>> img = nf.ni1.Nifti1Image(data, affine)
    >>> img.shape == (10, 1, 1)
    True
    >>> img2 = squeeze_image(img)
    >>> img2.shape == (10, 1, 1)
    True

    Only *final* dimensions of 1 are squeezed

    >>> shape = (1, 1, 5, 1, 2, 1, 1)
    >>> data = data.reshape(shape)
    >>> img = nf.ni1.Nifti1Image(data, affine)
    >>> img.shape == (1, 1, 5, 1, 2, 1, 1)
    True
    >>> img2 = squeeze_image(img)
    >>> img2.shape == (1, 1, 5, 1, 2)
    True
    """
    klass = img.__class__
    shape = img.shape
    slen = len(shape)
    if slen < 4:
        return klass.from_image(img)
    for bdim in shape[3::][::-1]:
        if bdim == 1:
            slen -= 1
        else:
            break
    if slen == len(shape):
        return klass.from_image(img)
    shape = shape[:slen]
    data = np.asanyarray(img.dataobj).reshape(shape)
    return klass(data,
                 img.affine,
                 img.header,
                 img.extra)

def concat_images(images, check_affines=True, axis=None):
    r""" Concatenate images in list to single image, along specified dimension

    Parameters
    ----------
    images : sequence
       sequence of ``SpatialImage`` or filenames of the same dimensionality\s
    check_affines : {True, False}, optional
       If True, then check that all the affines for `images` are nearly
       the same, raising a ``ValueError`` otherwise.  Default is True
    axis : None or int, optional
        If None, concatenates on a new dimension.  This requires all images to
        be the same shape.  If not None, concatenates on the specified
        dimension.  This requires all images to be the same shape, except on
        the specified dimension.

    Returns
    -------
    concat_img : ``SpatialImage``
       New image resulting from concatenating `images` across last
       dimension
    """
    images = [load(img) if not hasattr(img, 'get_data')
              else img for img in images]
    n_imgs = len(images)
    if n_imgs == 0:
        raise ValueError("Cannot concatenate an empty list of images.")
    img0 = images[0]
    affine = img0.affine
    header = img0.header
    klass = img0.__class__
    shape0 = img0.shape
    n_dim = len(shape0)
    if axis is None:
        # collect images in output array for efficiency
        out_shape = (n_imgs, ) + shape0
        out_data = np.empty(out_shape)
    else:
        # collect images in list for use with np.concatenate
        out_data = [None] * n_imgs
    # Get part of shape we need to check inside loop
    idx_mask = np.ones((n_dim,), dtype=bool)
    if axis is not None:
        idx_mask[axis] = False
    masked_shape = np.array(shape0)[idx_mask]
    for i, img in enumerate(images):
        if len(img.shape) != n_dim:
            raise ValueError(f'Image {i} has {len(img.shape)} dimensions, image 0 has {n_dim}')
        if not np.all(np.array(img.shape)[idx_mask] == masked_shape):
            raise ValueError(f'shape {img.shape} for image {i} not compatible with '
                             f'first image shape {shape0} with axis == {axis}')
        if check_affines and not np.all(img.affine == affine):
            raise ValueError(f'Affine for image {i} does not match affine for first image')
        # Do not fill cache in image if it is empty
        out_data[i] = np.asanyarray(img.dataobj)

    if axis is None:
        out_data = np.rollaxis(out_data, 0, out_data.ndim)
    else:
        out_data = np.concatenate(out_data, axis=axis)

    return klass(out_data, affine, header)

def four_to_three(img):
    """ Create 3D images from 4D image by slicing over last axis

    Parameters
    ----------
    img :  image
       4D image instance of some class with methods ``get_data``,
       ``header`` and ``affine``, and a class constructor
       allowing klass(data, affine, header)

    Returns
    -------
    imgs : list
       list of 3D images
    """
    arr = np.asanyarray(img.dataobj)
    header = img.header
    affine = img.affine
    image_maker = img.__class__
    if arr.ndim != 4:
        raise ValueError('Expecting four dimensions')
    imgs = []
    for i in range(arr.shape[3]):
        arr3d = arr[..., i]
        img3d = image_maker(arr3d, affine, header)
        imgs.append(img3d)
    return imgs

def as_closest_canonical(img, enforce_diag=False):
    """ Return `img` with data reordered to be closest to canonical

    Canonical order is the ordering of the output axes.

    Parameters
    ----------
    img : ``spatialimage``
    enforce_diag : {False, True}, optional
       If True, before transforming image, check if the resulting image
       affine will be close to diagonal, and if not, raise an error

    Returns
    -------
    canonical_img : ``spatialimage``
       Version of `img` where the underlying array may have been
       reordered and / or flipped so that axes 0,1,2 are those axes in
       the input data that are, respectively, closest to the output axis
       orientation.  We modify the affine accordingly.  If `img` is
       already has the correct data ordering, we just return `img`
       unmodified.
    """
    # Get the image class to transform the data for us
    img = img.as_reoriented(io_orientation(img.affine))

    # however, the affine may not be diagonal
    if enforce_diag and not _aff_is_diag(img.affine):
        raise OrientationError('Transformed affine is not diagonal')

    return img

def _aff_is_diag(aff):
    """ Utility function returning True if affine is nearly diagonal """
    rzs_aff = aff[:3, :3]
    return np.allclose(rzs_aff, np.diag(np.diag(rzs_aff)))


text:
# Copyright (c) 2020 PaddlePaddle Authors. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

from __future__ import print_function

import unittest
import paddle
import numpy as np
import paddle.fluid.core as core
from op_test import OpTest
import paddle.fluid as fluid
from paddle.fluid import Program, program_guard

class TestRepeatInterleaveOp(OpTest):

    def setUp(self):
        self.op_type = "repeat_interleave"
        self.init_dtype_type()
        index_np = np.random.randint(
            low=0, high=3, size=self.index_size).astype(self.index_type)
        x_np = np.random.random(self.x_shape).astype(self.x_type)
        self.inputs = {'X': x_np, 'RepeatsTensor': index_np}
        self.attrs = {'dim': self.dim}

        outer_loop = np.prod(self.x_shape[:self.dim])
        x_reshape = [outer_loop] + list(self.x_shape[self.dim:])
        x_np_reshape = np.reshape(x_np, tuple(x_reshape))
        out_list = []
        for i in range(outer_loop):
            for j in range(self.index_size):
                for k in range(index_np[j]):
                    out_list.append(x_np_reshape[i, j])
        self.out_shape = list(self.x_shape)
        self.out_shape[self.dim] = np.sum(index_np)
        self.out_shape = tuple(self.out_shape)

        out = np.reshape(out_list, self.out_shape)
        self.outputs = {'Out': out}

    def init_dtype_type(self):
        self.dim = 1
        self.x_type = np.float64
        self.index_type = np.int64
        self.x_shape = (8, 4, 5)
        self.index_size = self.x_shape[self.dim]

    def test_check_output(self):
        self.check_output()

    def test_check_grad_normal(self):
        self.check_grad(['X'], 'Out')

class TestRepeatInterleaveOp2(OpTest):

    def setUp(self):
        self.op_type = "repeat_interleave"
        self.init_dtype_type()
        index_np = 2
        x_np = np.random.random(self.x_shape).astype(self.x_type)
        self.inputs = {'X': x_np}  #, 'RepeatsTensor': None}
        self.attrs = {'dim': self.dim, 'Repeats': index_np}

        outer_loop = np.prod(self.x_shape[:self.dim])
        x_reshape = [outer_loop] + list(self.x_shape[self.dim:])
        x_np_reshape = np.reshape(x_np, tuple(x_reshape))
        out_list = []
        for i in range(outer_loop):
            for j in range(self.index_size):
                for k in range(index_np):
                    out_list.append(x_np_reshape[i, j])
        self.out_shape = list(self.x_shape)
        self.out_shape[self.dim] = index_np * self.index_size
        self.out_shape = tuple(self.out_shape)

        out = np.reshape(out_list, self.out_shape)
        self.outputs = {'Out': out}

    def init_dtype_type(self):
        self.dim = 1
        self.x_type = np.float64
        self.x_shape = (8, 4, 5)
        self.index_size = self.x_shape[self.dim]

    def test_check_output(self):
        self.check_output()

    def test_check_grad_normal(self):
        self.check_grad(['X'], 'Out')

class TestIndexSelectAPI(unittest.TestCase):

    def input_data(self):
        self.data_x = np.array([[1.0, 2.0, 3.0, 4.0], [5.0, 6.0, 7.0, 8.0],
                                [9.0, 10.0, 11.0, 12.0]])
        self.data_index = np.array([0, 1, 2, 1]).astype('int32')

    def test_repeat_interleave_api(self):
        paddle.enable_static()
        self.input_data()

        # case 1:
        with program_guard(Program(), Program()):
            x = fluid.layers.data(name='x', shape=[-1, 4])
            index = fluid.layers.data(name='repeats',
                                      shape=[4],
                                      dtype='int32',
                                      append_batch_size=False)
            z = paddle.repeat_interleave(x, index, axis=1)
            exe = fluid.Executor(fluid.CPUPlace())
            res, = exe.run(feed={
                'x': self.data_x,
                'repeats': self.data_index
            },
                           fetch_list=[z.name],
                           return_numpy=False)
        expect_out = np.repeat(self.data_x, self.data_index, axis=1)
        self.assertTrue(np.allclose(expect_out, np.array(res)))

        # case 2:
        repeats = np.array([1, 2, 1]).astype('int32')
        with program_guard(Program(), Program()):
            x = fluid.layers.data(name='x', shape=[-1, 4])
            index = fluid.layers.data(name='repeats',
                                      shape=[3],
                                      dtype='int32',
                                      append_batch_size=False)
            z = paddle.repeat_interleave(x, index, axis=0)
            exe = fluid.Executor(fluid.CPUPlace())
            res, = exe.run(feed={
                'x': self.data_x,
                'repeats': repeats,
            },
                           fetch_list=[z.name],
                           return_numpy=False)
        expect_out = np.repeat(self.data_x, repeats, axis=0)
        self.assertTrue(np.allclose(expect_out, np.array(res)))

        repeats = 2
        with program_guard(Program(), Program()):
            x = fluid.layers.data(name='x', shape=[-1, 4])
            z = paddle.repeat_interleave(x, repeats, axis=0)
            exe = fluid.Executor(fluid.CPUPlace())
            res, = exe.run(feed={'x': self.data_x},
                           fetch_list=[z.name],
                           return_numpy=False)
        expect_out = np.repeat(self.data_x, repeats, axis=0)
        self.assertTrue(np.allclose(expect_out, np.array(res)))

    def test_dygraph_api(self):
        self.input_data()
        # case axis none
        input_x = np.array([[1, 2, 1], [1, 2, 3]]).astype('int32')
        index_x = np.array([1, 1, 2, 1, 2, 2]).astype('int32')

        with fluid.dygraph.guard():
            x = fluid.dygraph.to_variable(input_x)
            index = fluid.dygraph.to_variable(index_x)
            z = paddle.repeat_interleave(x, index, None)
            np_z = z.numpy()
        expect_out = np.repeat(input_x, index_x, axis=None)
        self.assertTrue(np.allclose(expect_out, np_z))

        # case repeats int
        with fluid.dygraph.guard():
            x = fluid.dygraph.to_variable(input_x)
            index = 2
            z = paddle.repeat_interleave(x, index, None)
            np_z = z.numpy()
        expect_out = np.repeat(input_x, index, axis=None)
        self.assertTrue(np.allclose(expect_out, np_z))

        # case 1:
        with fluid.dygraph.guard():
            x = fluid.dygraph.to_variable(self.data_x)
            index = fluid.dygraph.to_variable(self.data_index)
            z = paddle.repeat_interleave(x, index, -1)
            np_z = z.numpy()
        expect_out = np.repeat(self.data_x, self.data_index, axis=-1)
        self.assertTrue(np.allclose(expect_out, np_z))

        with fluid.dygraph.guard():
            x = fluid.dygraph.to_variable(self.data_x)
            index = fluid.dygraph.to_variable(self.data_index)
            z = paddle.repeat_interleave(x, index, 1)
            np_z = z.numpy()
        expect_out = np.repeat(self.data_x, self.data_index, axis=1)
        self.assertTrue(np.allclose(expect_out, np_z))

        # case 2:
        index_x = np.array([1, 2, 1]).astype('int32')
        with fluid.dygraph.guard():
            x = fluid.dygraph.to_variable(self.data_x)
            index = fluid.dygraph.to_variable(index_x)
            z = paddle.repeat_interleave(x, index, axis=0)
            np_z = z.numpy()
        expect_out = np.repeat(self.data_x, index, axis=0)
        self.assertTrue(np.allclose(expect_out, np_z))

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


text:
import cv2
import glob
import pickle
import numpy as np
import matplotlib.image as mpimg
from skimage.feature import hog
from scipy.ndimage.measurements import label
#from os import walk
from os import path
import time
from sklearn.svm import LinearSVC
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from scipy.ndimage.measurements import label
from copy import copy
#from timeit import default_timer as timer

### Parameters
color_space = 'YCrCb' # Can be RGB, HSV, LUV, HLS, YUV, YCrCb
orient = 9  # HOG orientations
pix_per_cell = 8 # 16 # 8 # HOG pixels per cell
cell_per_block = 1 # 2 # HOG cells per block, which can handel e.g. shadows
hog_channel = 0 #'ALL' # Can be 0, 1, 2, or 'ALL'
spatial_size = (16, 16) # (32, 32) # Spatial binning dimensions
hist_bins = 32    # Number of histogram bins
spatial_feat = True # True # Spatial features on or off
hist_feat = True # True # Histogram features on or off
hog_feat = True # HOG features on or off
clf_path = 'car_svc_pickle.p'  # if classifier exist

Heatmap_buffer = []
N_buffer = 10 #3
buffer_weights=np.arange(1,N_buffer+1)/N_buffer

xy_windows = [(128, 128), (96, 96), (80, 80)]
xy_overlaps= [(0.5, 0.5), (0.5, 0.5), (0.5, 0.5)]
x_start_stops = [[None, None], [32, None], [400, 1280]]
y_start_stops = [[400, 650], [400, 650], [390, 550]]

# y_start_stop = [400, 656] # Min and max in y to search in slide_window()
# ystart_0 = y_start_stop[0]
# ystop_0 = ystart_0 + 64*2
# ystart_1 = ystart_0
# ystop_1 = y_start_stop[1]
# ystart_2 = ystart_0
# ystop_2 = y_start_stop[1]
# ystarts = [ystart_1, ystart_2]
# ystops = [ystop_1-100, ystop_2]
# search_window_scales = [1.5, 2]  # (64x64), (96x96), (128x128)

class AverageHotBox ():
    """Class that covers joining hot boxes algorithm.
    Idea is take fist box (called average box) form input boxes
    and join it with all boxes that are close enough.
    After joining two boxes we need to update average box (here
    just increasing size to cover both joining boxes).
    Loop until cannot further join.
    For left boxes repeat all procedure.
    As a result we also get average boxes strengths - the number of
    boxes it was joined to.
    This class represents one average box.
    """
    def __init__ (self, box):
        self.avg_box = [list(p) for p in box]
        self.detected_count = 1
        self.boxes = [box]

    def get_strength (self):
        """Returns number of joined boxes"""
        return self.detected_count

    def get_box (self):
        """Uses joined boxes information to compute
        this average box representation as hot box.
        This box has average center of all boxes and have
        size of 2 standard deviation by x and y coordinates of its points
        """
        if len(self.boxes) > 1:
            center = np.average (np.average (self.boxes, axis=1), axis=0).astype(np.int32).tolist()

            # getting all x and y coordinates of
            # all corners of joined boxes separately
            xs = np.array(self.boxes) [:,:,0]
            ys = np.array(self.boxes) [:,:,1]

            half_width = int(np.std (xs))
            half_height = int(np.std (ys))
            return (
                (
                    center[0] - half_width,
                    center[1] - half_height
                ), (
                    center[0] + half_width,
                    center[1] + half_height
                ))
        else:
            return self.boxes [0]

    def is_close (self, box):
        """Check wether specified box is close enough for joining
        to be close need to overlap by 30% of area of this box or the average box
        """
        # Reference -
        # http://math.stackexchange.com/questions/99565/simplest-way-to-calculate-the-intersect-area-of-two-rectangles
        x11 = self.avg_box [0][0]
        y11 = self.avg_box [0][1]
        x12 = self.avg_box [1][0]
        y12 = self.avg_box [1][1]
        x21 = box [0][0]
        y21 = box [0][1]
        x22 = box [1][0]
        y22 = box [1][1]

        x_overlap = max(0, min(x12,x22) - max(x11,x21))
        y_overlap = max(0, min(y12,y22) - max(y11,y21))

        area1 = (x12 - x11) * (y12 - y11)
        area2 = (x22 - x21) * (y22 - y21)
        intersection = x_overlap * y_overlap;

        if (
            intersection >= 0.3 * area1 or
            intersection >= 0.3 * area2
        ):
            return True
        else:
            return False

    def join (self, boxes):
        """Join in all boxes from list of given boxes,
        removes joined boxes from input list of boxes
        """

        joined = False

        for b in boxes:
            if self.is_close (b):
                boxes.remove (b)
                self.boxes.append (b)
                self.detected_count += 1

                self.avg_box [0][0] = min (self.avg_box [0][0], b [0][0])
                self.avg_box [0][1] = min (self.avg_box [0][1], b [0][1])
                self.avg_box [1][0] = max (self.avg_box [1][0], b [1][0])
                self.avg_box [1][1] = max (self.avg_box [1][1], b [1][1])

                joined = True

        return joined

def calc_average_boxes (hot_boxes, strength):
    """Compute average boxes from specified hot boxes and returns
    average boxes with equals or higher strength
    """
    avg_boxes = []
    while len(hot_boxes) > 0:
        b = hot_boxes.pop (0)
        hb = AverageHotBox (b)
        while hb.join (hot_boxes):
            pass
        avg_boxes.append (hb)

    boxes = []
    for ab in avg_boxes:
        if ab.get_strength () >= strength:
            boxes.append (ab.get_box ())
    return boxes

class LastHotBoxesQueue ():
    """Class for accumulation of hot boxes from last 10 frames
    """
    def __init__ (self):
        self.queue_max_len = 10 # number items to store
        self.last_boxes = []

    def put_hot_boxes (self, boxes):
        """Put frame hot boxes
        """
        if (len(self.last_boxes) > self.queue_max_len):
            tmp = self.last_boxes.pop (0)

        self.last_boxes.append (boxes)

    def get_hot_boxes (self):
        """Get last 10 frames hot boxes
        """
        b = []
        for boxes in self.last_boxes:
            b.extend (boxes)
        return b

last_hot_boxes = LastHotBoxesQueue()
# heat_map = deque(np.array([np.zeros(img_size).astype(np.float)]), maxlen=N_buffer)
def car_svm(image_in): #, img_lane_augmented, lane_info):

    #start = timer()
    image = np.copy(image_in)
    img = image.astype(np.float32)/255

    windows_list, img_box = hot_windows(img)

    last_hot_boxes.put_hot_boxes(windows_list)
    hot_boxes = last_hot_boxes.get_hot_boxes ()

    # calculating average boxes and use strong ones
    # need to tune strength on particular classifer
    avg_boxes = calc_average_boxes (hot_boxes, 12) # 20)
    # heatmap = generate_heatmap(img, avg_boxes)
    # labels = label(heatmap)
    # cv2.putText(annotated_image1,"Close Vehicles: " + "{:0.0f}".format(labels[1]), org=(50,200), fontFace=cv2.FONT_HERSHEY_SIMPLEX,
    #     fontScale=1, color=(255,255,255), lineType = cv2.LINE_AA, thickness=2)
    # image_with_boxes = draw_labeled_bboxes(img, labels)

    image_with_boxes = draw_boxes(img, avg_boxes, color=(0, 0, 1), thick=4)
    bboxes=[]
    for ab in avg_boxes:
        ((x1,y1),(x2,y2)) = ab
        bboxes.append([x1,y1,x2,y2])
    return image_with_boxes * 255, bboxes

    # for search_window_scale, ystart, ystop in zip(search_window_scales, ystarts, ystops):
    #     windows_list_tmp = find_cars(np.copy(image), ystart, ystop, search_window_scale, svc, X_scaler, orient, pix_per_cell, cell_per_block,
    #                         spatial_size, hist_bins, hog_channel, color_space, spatial_feat, hist_feat, hog_feat)
    #     windows_list.extend(windows_list_tmp)

    # heatmap_pre = generate_heatmap(img, windows_list)
    # draw_img, heatmap_post, bboxes = draw_bboxes(image, copy(Heatmap_buffer), heatmap_pre, min(len(Heatmap_buffer)+1,N_buffer) )

    # if len(Heatmap_buffer) >= N_buffer:
    #     Heatmap_buffer.pop(0)

    # return draw_img

def draw_labeled_bboxes(img, labels):
    # Iterate through all detected cars
    for car_number in range(1, labels[1]+1):
        # Find pixels with each car number label value
        nonzero = (labels[0] == car_number).nonzero()
        # Identify x and y values of those pixels
        nonzeroy = np.array(nonzero[0])
        nonzerox = np.array(nonzero[1])
        # Define a bounding box based on min/max x and y
        min_x = np.min(nonzerox)
        min_y = np.min(nonzeroy)
        max_x = np.max(nonzerox)
        max_y = np.max(nonzeroy)
        if (max_x - min_x) > 30 or (max_y - max_y)> 30:
            bbox = ((np.min(nonzerox), np.min(nonzeroy)), (np.max(nonzerox), np.max(nonzeroy)))
            # Draw the box on the image
            cv2.rectangle(img, bbox[0], bbox[1],(0,0,255), 6)
            cv2.putText(img,"Vehicle: " + "{:0.0f}".format(car_number), org=(bbox[0][0],bbox[0][1]-20), fontFace=cv2.FONT_HERSHEY_SIMPLEX,
                fontScale=0.7, color=(255,255,255), lineType = cv2.LINE_AA, thickness=2)
    # return the image
    return img

# Define the hot windows by slide window search
# Returns: hot_windows(list), image_with_hot_windows_drawn(numpy.array)
def hot_windows(image):
    dst = np.copy (image)
    all_hot_wins = []

    # iterate over previousely defined sliding windows
    for x_limits, y_limits, window_size, overlap in zip(x_start_stops, y_start_stops, xy_windows, xy_overlaps):
        windows = slide_window( dst,
            x_start_stop=x_limits,
            y_start_stop=y_limits,
            xy_window=window_size,
            xy_overlap=overlap)

        hot_wins = search_windows(image, windows, svc, X_scaler, color_space=color_space,
                    spatial_size=spatial_size, hist_bins=hist_bins,
                    orient=orient, pix_per_cell=pix_per_cell,
                    cell_per_block=cell_per_block,
                    hog_channel=hog_channel, spatial_feat=spatial_feat,
                    hist_feat=hist_feat, hog_feat=hog_feat)

        all_hot_wins.extend(hot_wins)
        dst = draw_boxes(dst, hot_wins, color=(0, 0, 1), thick=4)

    return all_hot_wins, dst

# Define a function to return HOG features and visualization
def get_hog_features(img, orient, pix_per_cell, cell_per_block,
                        vis=False, feature_vec=True):
    # Call with two outputs if vis==True
    if vis == True:
        features, hog_image = hog(img, orientations=orient,
                                  pixels_per_cell=(pix_per_cell, pix_per_cell),
                                  cells_per_block=(cell_per_block, cell_per_block),
                                  transform_sqrt=True,
                                  visualise=vis, feature_vector=feature_vec)
        return features, hog_image
    # Otherwise call with one output
    else:
        features = hog(img, orientations=orient,
                       pixels_per_cell=(pix_per_cell, pix_per_cell),
                       cells_per_block=(cell_per_block, cell_per_block),
                       transform_sqrt=True,
                       visualise=vis, feature_vector=feature_vec)
        return features

def bin_spatial(img, size=(32, 32)):
    color1 = cv2.resize(img[:,:,0], size).ravel()
    color2 = cv2.resize(img[:,:,1], size).ravel()
    color3 = cv2.resize(img[:,:,2], size).ravel()
    return np.hstack((color1, color2, color3))

# Define a function to compute color histogram features
# NEED TO CHANGE bins_range if reading .png files with mpimg!
def color_hist(img, nbins=32):
    # Compute the histogram of the color channels separately
    channel1_hist = np.histogram(img[:,:,0], bins=nbins)
    channel2_hist = np.histogram(img[:,:,1], bins=nbins)
    channel3_hist = np.histogram(img[:,:,2], bins=nbins)
    # Concatenate the histograms into a single feature vector
    hist_features = np.concatenate((channel1_hist[0], channel2_hist[0], channel3_hist[0]))
    # Return the individual histograms, bin_centers and feature vector
    return hist_features

# convert from RGB to another color space
def convert_color(img, color='YCrCb'):
    if color != 'RGB':
        image = cv2.cvtColor(img, eval('cv2.COLOR_RGB2'+color))
    else: image = np.copy(img)
    return image

# Define a function to extract features from a list of images
# Have this function call bin_spatial() and color_hist()
def extract_features(imgs, color_space='RGB', spatial_size=(32, 32),
                        hist_bins=32, orient=9,
                        pix_per_cell=8, cell_per_block=2, hog_channel=0,
                        spatial_feat=True, hist_feat=True, hog_feat=True):
    # Create a list to append feature vectors to
    features = []
    # Iterate through the list of images
    for file in imgs:
        file_features = []

        # png is scale from (0,1)
        image = mpimg.imread(file)

        # apply color conversion if other than 'RGB'
        feature_image = convert_color(image, color_space)

        if spatial_feat == True:
            spatial_features = bin_spatial(feature_image, size=spatial_size)
            file_features.append(spatial_features)
        if hist_feat == True:
            # Apply color_hist()
            hist_features = color_hist(feature_image, nbins=hist_bins)
            file_features.append(hist_features)
        if hog_feat == True:
        # Call get_hog_features() with vis=False, feature_vec=True
            if hog_channel == 'ALL':
                hog_features = []
                for channel in range(feature_image.shape[2]):
                    hog_features.append(get_hog_features(feature_image[:,:,channel],
                                        orient, pix_per_cell, cell_per_block,
                                        vis=False, feature_vec=True))
                hog_features = np.ravel(hog_features)
            else:
                hog_features = get_hog_features(feature_image[:,:,hog_channel], orient,
                            pix_per_cell, cell_per_block, vis=False, feature_vec=True)
            # Append the new feature vector to the features list
            file_features.append(hog_features)
        features.append(np.concatenate(file_features))
    # Return list of feature vectors
    return features

def slide_window(img, x_start_stop=[None, None], y_start_stop=[None, None],
                    xy_window=(64, 64), xy_overlap=(0.5, 0.5)):
    # If x and/or y start/stop positions not defined, set to image size
    if x_start_stop[0] == None:
        x_start_stop[0] = 0
    if x_start_stop[1] == None:
        x_start_stop[1] = img.shape[1]
    if y_start_stop[0] == None:
        y_start_stop[0] = 0
    if y_start_stop[1] == None:
        y_start_stop[1] = img.shape[0]
    # Compute the span of the region to be searched
    xspan = x_start_stop[1] - x_start_stop[0]
    yspan = y_start_stop[1] - y_start_stop[0]
    # Compute the number of pixels per step in x/y
    nx_pix_per_step = np.int(xy_window[0]*(1 - xy_overlap[0]))
    ny_pix_per_step = np.int(xy_window[1]*(1 - xy_overlap[1]))
    # Compute the number of windows in x/y
    nx_buffer = np.int(xy_window[0]*(xy_overlap[0]))
    ny_buffer = np.int(xy_window[1]*(xy_overlap[1]))
    nx_windows = np.int((xspan-nx_buffer)/nx_pix_per_step)
    ny_windows = np.int((yspan-ny_buffer)/ny_pix_per_step)
    # Initialize a list to append window positions to
    window_list = []
    for ys in range(ny_windows):
        for xs in range(nx_windows):
            # Calculate window position
            startx = xs*nx_pix_per_step + x_start_stop[0]
            endx = startx + xy_window[0]
            starty = ys*ny_pix_per_step + y_start_stop[0]
            endy = starty + xy_window[1]

            # Append window position to list
            window_list.append(((startx, starty), (endx, endy)))
    # Return the list of windows
    return window_list

# Define a function to draw bounding boxes
def draw_boxes(img, bboxes, color=(0, 0, 255), thick=6):
    # Make a copy of the image
    imcopy = np.copy(img)
    # Iterate through the bounding boxes
    for bbox in bboxes:
        # Draw a rectangle given bbox coordinates
        cv2.rectangle(imcopy, bbox[0], bbox[1], color, thick)
    # Return the image copy with boxes drawn
    return imcopy

# Define a function to extract features from a single image window
def single_img_features(img, color_space='RGB', spatial_size=(32, 32),
                        hist_bins=32, orient=9,
                        pix_per_cell=8, cell_per_block=2, hog_channel=0,
                        spatial_feat=True, hist_feat=True, hog_feat=True):
    #1) Define an empty list to receive features
    img_features = []
    #2) Apply color conversion if other than 'RGB'
    feature_image = convert_color(img, color_space)
    #3) Compute spatial features if flag is set
    if spatial_feat == True:
        spatial_features = bin_spatial(feature_image, size=spatial_size)
        #4) Append features to list
        img_features.append(spatial_features)
    #5) Compute histogram features if flag is set
    if hist_feat == True:
        hist_features = color_hist(feature_image, nbins=hist_bins)
        #6) Append features to list
        img_features.append(hist_features)
    #7) Compute HOG features if flag is set
    if hog_feat == True:
        if hog_channel == 'ALL':
            hog_features = []
            for channel in range(feature_image.shape[2]):
                hog_features.extend(get_hog_features(feature_image[:,:,channel],
                                    orient, pix_per_cell, cell_per_block,
                                    vis=False, feature_vec=True))
        else:
            hog_features = get_hog_features(feature_image[:,:,hog_channel], orient,
                        pix_per_cell, cell_per_block, vis=False, feature_vec=True)
        #8) Append features to list
        img_features.append(hog_features)

    #9) Return concatenated array of features
    return np.concatenate(img_features)

# Define a function you will pass an image
# and the list of windows to be searched (output of slide_windows())
def search_windows(img, windows, clf, scaler, color_space='RGB',
                    spatial_size=(32, 32), hist_bins=32, orient=9,
                    pix_per_cell=8, cell_per_block=2,
                    hog_channel=0, spatial_feat=True,
                    hist_feat=True, hog_feat=True):

    #1) Create an empty list to receive positive detection windows
    on_windows = []
    #2) Iterate over all windows in the list
    for window in windows:
        #3) Extract the test window from original image
        test_img = cv2.resize(img[window[0][1]:window[1][1], window[0][0]:window[1][0]], (64, 64))  # training image is (64,64)

        #4) Extract features for that window using single_img_features()
        features = single_img_features(test_img, color_space=color_space,
                            spatial_size=spatial_size, hist_bins=hist_bins,
                            orient=orient, pix_per_cell=pix_per_cell,
                            cell_per_block=cell_per_block,
                            hog_channel=hog_channel, spatial_feat=spatial_feat,
                            hist_feat=hist_feat, hog_feat=hog_feat)
        #5) Scale extracted features to be fed to classifier
        X = np.array(features).reshape(1, -1)
        test_features = scaler.transform(X)
        #6) Predict using your classifier
        prediction = clf.predict(test_features)
        #7) If positive (prediction == 1) then save the window
        if prediction == 1:
            on_windows.append(window)
    #8) Return windows for positive detections
    return on_windows

def find_cars(img, ystart, ystop, scale, svc, X_scaler, orient, pix_per_cell, cell_per_block, spatial_size, hist_bins,
    hog_channel, color_space, spatial_feat, hist_feat, hog_feat): # , dec_func_thresh=1.0):

    on_windows = []
    features = []
    img = img.astype(np.float32)/255

    img_tosearch = img[ystart:ystop,:,:]
    ctrans_tosearch = convert_color(img_tosearch, color=color_space)

    if scale != 1:
        imshape = ctrans_tosearch.shape
        ctrans_tosearch = cv2.resize(ctrans_tosearch, (np.int(imshape[1]/scale), np.int(imshape[0]/scale)))

    ch1 = ctrans_tosearch[:,:,0]
    ch2 = ctrans_tosearch[:,:,1]
    ch3 = ctrans_tosearch[:,:,2]

    # Define blocks and steps as above, hold the number of hog cells
    nxblocks = (ch1.shape[1] // pix_per_cell) - cell_per_block + 1
    nyblocks = (ch1.shape[0] // pix_per_cell) - cell_per_block + 1
    nfeat_per_block = orient*cell_per_block**2

    # 64 was the orginal sampling rate, with 8 cells and 8 pix per cell
    window = 64
    nblocks_per_window = (window // pix_per_cell) - cell_per_block + 1
    cells_per_step = 2  # Instead of overlap, define how many cells to step: there are 8 cells, and move 2 cells per step, 75% overlap
    nxsteps = (nxblocks - nblocks_per_window) // cells_per_step + 1
    nysteps = (nyblocks - nblocks_per_window) // cells_per_step + 1

    # Compute individual channel HOG features for the entire image
    hog1 = get_hog_features(ch1, orient, pix_per_cell, cell_per_block, feature_vec=False)
    hog2 = get_hog_features(ch2, orient, pix_per_cell, cell_per_block, feature_vec=False)
    hog3 = get_hog_features(ch3, orient, pix_per_cell, cell_per_block, feature_vec=False)

    for xb in range(nxsteps):
        for yb in range(nysteps):
            ypos = yb*cells_per_step
            xpos = xb*cells_per_step

            if hog_feat:
                # Extract HOG for this patch
                if hog_channel == 0:
                    hog_features = hog1[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                elif hog_channel == 1:
                    hog_features = hog2[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                elif hog_channel == 2:
                    hog_features = hog3[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                else:
                    hog_feat1 = hog1[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                    hog_feat2 = hog2[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                    hog_feat3 = hog3[ypos:ypos+nblocks_per_window, xpos:xpos+nblocks_per_window].ravel()
                    hog_features = np.hstack((hog_feat1, hog_feat2, hog_feat3))

            xleft = xpos*pix_per_cell
            ytop = ypos*pix_per_cell

            features.append(hog_features)
            # Extract the image patch
            subimg = cv2.resize(ctrans_tosearch[ytop:ytop+window, xleft:xleft+window], (64,64))

            # Get color features
            if spatial_feat:
                spatial_features = bin_spatial(subimg, size=spatial_size)
                features.append(spatial_feat)

            if hist_feat:
                hist_features = color_hist(subimg, nbins=hist_bins)
                features.append(hist_features)

            #X = np.hstack((spatial_features, hist_features, hog_features)).reshape(1, -1)
            X = np.array(np.concatenate(features)).reshape(1, -1)

            test_features = X_scaler.transform(X)
            test_prediction = svc.predict(test_features)
            # test_dec_func = svc.decision_function(test_features)

            if test_prediction == 1 : # and test_dec_func >= dec_func_thresh :
                xbox_left = np.int(xleft*scale)
                ytop_draw = np.int(ytop*scale)
                win_draw = np.int(window*scale)
                # cv2.rectangle(draw_img,(xbox_left, ytop_draw+ystart),(xbox_left+win_draw,ytop_draw+win_draw+ystart),(0,0,255),6)
                on_windows.append(((xbox_left, ytop_draw+ystart),(xbox_left+win_draw,ytop_draw+win_draw+ystart)))

    return on_windows

def add_heat(heatmap, bbox_list):
    # Iterate through list of bboxes
    for box in bbox_list:
        # Add += 1 for all pixels inside each bbox
        # Assuming each "box" takes the form ((x1, y1), (x2, y2))
        heatmap[box[0][1]:box[1][1], box[0][0]:box[1][0]] += 1

    # Return updated heatmap
    return heatmap

def apply_threshold(heatmap, threshold):
    # Zero out pixels below the threshold
    heatmap[heatmap <= threshold] = 0
    # Return thresholded map
    return heatmap

def draw_bboxes(img, heatmap_buffer, heatmap_pre, N_buffer):

    heatmap_buffer.append(heatmap_pre)

    if len(heatmap_buffer) > N_buffer: # remove the first component if it is more than N_buffer elements
        heatmap_buffer.pop(0)

    # weight the heatmap based on current frame and previous N frames
    idxs = range(N_buffer)
    for b, w, idx in zip(heatmap_buffer, buffer_weights, idxs):
        heatmap_buffer[idx] = b * w

    heatmap = np.sum(np.array(heatmap_buffer), axis=0)
    heatmap = apply_threshold( heatmap, threshold= sum(buffer_weights[0:N_buffer])*2)

    # Find final boxes from heatmap using label function
    labels = label(heatmap)

    bboxes = []
    # locate the bounding box
    for car_number in range(1, labels[1]+1):
        # Find pixels with each car_number label value
        nonzero = (labels[0] == car_number).nonzero()
        # Identify x and y values of those pixels
        nonzeroy = np.array(nonzero[0])
        nonzerox = np.array(nonzero[1])
        # Define a bounding box based on min/max x and y
        bbox_tmp = ((np.min(nonzerox), np.min(nonzeroy)), (np.max(nonzerox), np.max(nonzeroy)))
        bboxes.append(bbox_tmp)

    for bbox in bboxes:
        # Draw the box on the image
        cv2.rectangle(img, bbox[0], bbox[1], (0,0,255), 4)

    # Return the image

    return img, heatmap, bboxes

def generate_heatmap(image, windows_list):
    heat = np.zeros_like(image[:,:,0]).astype(np.float)

    # Add heat to each box in box list
    heat = add_heat(heat, windows_list)

    # Apply threshold to help remove false positives
    heat = apply_threshold(heat, 1)

    # Visualize the heatmap when displaying
    heatmap = np.clip(heat, 0, 255)

    # if np.amax(heatmap) == 0:
    #     flag = False
    # else:
    #     flag = True

    return heatmap

# if svm classifer exist, load it; otherwise, compute the svm classifier
if path.isfile(clf_path):

    print('loading existing classifier...')
    with open(clf_path, 'rb') as file:
        clf_pickle = pickle.load(file)
        svc = clf_pickle["svc"]
        X_scaler = clf_pickle["scaler"]
        orient = clf_pickle["orient"]
        pix_per_cell = clf_pickle["pix_per_cell"]
        cell_per_block = clf_pickle["cell_per_block"]
        spatial_size = clf_pickle["spatial_size"]
        hist_bins = clf_pickle["hist_bins"]
        color_space = clf_pickle["color_space"]

        hog_channel = clf_pickle["hog_channel"]
        spatial_feat = clf_pickle["spatial_feat"]
        hist_feat = clf_pickle["hist_feat"]
        hog_feat = clf_pickle["hog_feat"]
else:
    # Read in cars and notcars
    cars = glob.glob('train_images/vehicles/**/*.png')
    notcars = glob.glob('train_images/non-vehicles/**/*.png')

    # set the sample size
    sample_size = min(len(cars), len(notcars))
    cars = cars[0:sample_size]
    notcars = notcars[0:sample_size]

    print('Parameter file', clf_path, 'does not exist, train the car model...')

    print('extracting car features...')
    car_features = extract_features(cars, color_space=color_space,
                            spatial_size=spatial_size, hist_bins=hist_bins,
                            orient=orient, pix_per_cell=pix_per_cell,
                            cell_per_block=cell_per_block,
                            hog_channel=hog_channel, spatial_feat=spatial_feat,
                            hist_feat=hist_feat, hog_feat=hog_feat)
    print('car features extracted!')
    print('extracting noncar features...')
    notcar_features = extract_features(notcars, color_space=color_space,
                            spatial_size=spatial_size, hist_bins=hist_bins,
                            orient=orient, pix_per_cell=pix_per_cell,
                            cell_per_block=cell_per_block,
                            hog_channel=hog_channel, spatial_feat=spatial_feat,
                            hist_feat=hist_feat, hog_feat=hog_feat)
    print('noncar features extracted!')

    X = np.vstack((car_features, notcar_features)).astype(np.float64)
    # Fit a per-column scaler
    X_scaler = StandardScaler().fit(X)
    # Apply the scaler to X
    scaled_X = X_scaler.transform(X)

    # Define the labels vector
    y = np.hstack((np.ones(len(car_features)), np.zeros(len(notcar_features))))

    # Split up data into randomized training and test sets
    # rand_state = np.random.randint(0, 100)
    X_train, X_test, y_train, y_test = train_test_split(
        scaled_X, y, test_size=0.2, random_state=40)

    print('Using:',orient,'orientations',pix_per_cell,
        'pixels per cell and', cell_per_block,'cells per block')
    print('Feature vector length:', len(X_train[0]))
    # Use a linear SVC
    svc =  LinearSVC(C=0.004) # , loss='hinge')
    # Check the training time for the SVC
    t=time.time()
    svc.fit(X_train, y_train)
    t2 = time.time()
    print(round(t2-t, 2), 'Seconds to train SVC...')
    # Check the score of the SVC
    print('Test Accuracy of SVC = ', round(svc.score(X_test, y_test), 4))
    # Check the prediction time for a single sample
    t=time.time()

    # save classifier
    clf_pickle = {}
    clf_pickle["svc"] = svc
    clf_pickle["scaler"] = X_scaler
    clf_pickle["orient"] = orient
    clf_pickle["pix_per_cell"] = pix_per_cell
    clf_pickle["cell_per_block"] = cell_per_block
    clf_pickle["spatial_size"] = spatial_size
    clf_pickle["hist_bins"] = hist_bins
    clf_pickle["color_space"] = color_space

    clf_pickle["hog_channel"] = hog_channel
    clf_pickle["spatial_feat"] = spatial_feat
    clf_pickle["hist_feat"] = hist_feat
    clf_pickle["hog_feat"] = hog_feat

    pickle.dump(clf_pickle, open(clf_path, "wb"))
    print("Classifier is written into: {}".format(clf_path))

print('svc = ',svc)
print('scaler = ',X_scaler)
print('orient = ',orient)
print('pix_per_cell = ',pix_per_cell)
print('cell_per_block = ',cell_per_block)
print('spatial_size = ',spatial_size)
print('hist_bins = ',hist_bins)
print('color_space = ',color_space)
print('hog_channel = ',hog_channel)
print('spatial_feat = ',spatial_feat)
print('hist_feat = ',hist_feat)
print('hog_feat = ',hog_feat)


text:
# -*- coding: utf-8 -*-
"""
@author: nicolas.posocco
"""

from .....utils import bounded_1d_kde
import numpy as np

def binary_reliability_curve(model=None, X=None, Y=None,
                             kernel=None, bandwidth=None,
                             positive_scores=None,
                             positive_scores_for_positive_gt=None,
                             positive_class_probability=None):
    if positive_scores is None:
        # model
        assert model is not None, "model has to be provided when positive_scores is not."
        assert hasattr(model, "predict_proba")

        # X
        assert X is not None, "X has to be provided when positive_scores is not."
        assert type(X) is np.ndarray
        assert X.ndim == 2

        positive_scores = model.predict_proba(X)[:, 1]

    if positive_class_probability is None:
        # Y
        assert Y is not None, "Y has to be provided when positive_class_probability is not."
        assert type(Y) is np.ndarray
        assert Y.ndim == 1

        # model
        assert model is not None, "model has to be provided when positive_scores is not."
        assert hasattr(model, "classes_")

        unique, counts = np.unique(Y, return_counts=True)
        counts = dict(zip(unique, counts))
        positive_class_probability = counts[model.classes_[1]] / Y.shape[0]

    if positive_scores_for_positive_gt is None:
        # Y
        assert Y is not None, "Y has to be provided when positive_scores_for_positive_gt is not."
        assert type(Y) is np.ndarray
        assert Y.ndim == 1

        # model
        assert model is not None, "model has to be provided when positive_scores_for_positive_gt is not."
        assert hasattr(model, "classes_")

        positive_scores_for_positive_gt = positive_scores[Y == model.classes_[1]]

    # kernel
    assert kernel is not None

    # Defining the grid on which the density is estimated
    # (its support is [0, 1], yet we need to perform mirroring to constraint the kde to this domain)
    grid = np.linspace(-3, 4, 2 ** 20)

    # Estimating the density of P( S+ )
    s, ps_y = bounded_1d_kde(positive_scores,
                             low_bound=0, high_bound=1,
                             kernel=kernel, bandwidth=bandwidth,
                             output="discrete_signal",
                             grid=grid)

    # Rectification to avoid division by zero bellow
    epsilon = 1e-8
    ps_y[ps_y < epsilon] = epsilon

    # Estimating the density of P( S+ | Y=+ )
    _, psky1_y = bounded_1d_kde(positive_scores_for_positive_gt,
                                low_bound=0, high_bound=1,
                                kernel=kernel, bandwidth=bandwidth,
                                output="discrete_signal",
                                grid=grid)

    # Estimating the density of P( Y=+ | S+ )
    # using bayes inversion of the density
    py1ks_y = psky1_y * positive_class_probability / ps_y

    # Rectification : zero score density -> no calibration error prior
    py1ks_y[ps_y == epsilon] = s[ps_y == epsilon]

    reliability_curve = {"scores": s, "errors": py1ks_y}

    return reliability_curve


text:
#!/usr/bin/python3

import os
os.environ['TF_CPP_MIN_LOG_LEVEL'] = '2'
import tensorflow as tf
tf.compat.v1.logging.set_verbosity(tf.compat.v1.logging.ERROR)

import model
import time
from PIL import Image, ImageFilter,ImageEnhance
import random
import numpy as np
import sys
from tqdm import tqdm
import math

IMG_SIZE = 1024
all_outputs = 1
pad_width = 24
total_steps = 0
progress_bar = None
current_progress= 0

args = sys.argv
if len(args) != 3:
  print("Wrong arguments. Usage  removestars.py input_file output_file")
  print ('Argument List:', str(sys.argv))
  exit(1)

def process_channel(channel,pad_width,input_image_size):
    global progress_bar,step_size,current_progress
    output_image =  Image.new('L', input_image_size)
    for i in range(0,int(channel.size[0]/IMG_SIZE)):
        for j in range(0,int(channel.size[1]/IMG_SIZE)):
            corp_rect = (i*IMG_SIZE,j*IMG_SIZE,i*IMG_SIZE+IMG_SIZE,j*IMG_SIZE+IMG_SIZE)
            current_tile = channel.crop(corp_rect)
            current_tile = current_tile.convert('L')
            blank_image =  current_tile.copy()
            current_tile = current_tile.resize((IMG_SIZE-pad_width*2,IMG_SIZE-pad_width*2))
            blank_image.paste(current_tile,(pad_width,pad_width))
            blank_image  = np.asarray(blank_image,dtype="float16").reshape(1,IMG_SIZE,IMG_SIZE)/255
            predicted_section = G2.predict(blank_image)
            predicted_section = predicted_section.reshape(IMG_SIZE,IMG_SIZE)*255
            predicted_section = Image.fromarray(predicted_section).convert('L')
            predicted_section = predicted_section.crop((pad_width,pad_width,IMG_SIZE-pad_width,IMG_SIZE-pad_width))
            predicted_section = predicted_section.resize((IMG_SIZE,IMG_SIZE))
            output_image.paste(predicted_section, (i*IMG_SIZE,j*IMG_SIZE),  mask=None)
            progress_bar.update(step_size)
            current_progress = current_progress + step_size
    return output_image

G2 = model.Generator()
G2.load_weights("weights/weights")

source_image = Image.open(args[1])
mode = source_image.mode
size = source_image.size
max_dimension = max(size)

a,b = divmod(max_dimension,1024)
if b > 0 :
  max_dimension = (a+1)*1024

input_image = source_image.crop((0,0,max_dimension,max_dimension))
progress_bar = tqdm(total=100)

if mode == 'L':
  total_steps = math.ceil((max_dimension/1024)**2)
  step_size = math.ceil((1/total_steps)*100)
  output = process_channel(input_image,pad_width,input_image.size)
  output = output.crop((0,0,size[0],size[1]))
  output.save(args[2])
  progress_bar.update(100-current_progress)
elif mode == 'RGB' or mode == 'RGBA':
  total_steps = int(3*(max_dimension/1024)**2)
  step_size = int((1/total_steps)*100)
  channels = Image.Image.split(input_image)
  all_outputs  = []
  for channel in channels[0:3]:
    channel_output = process_channel(channel,pad_width,input_image.size)
    channel_output = channel_output.crop((0,0,size[0],size[1]))
    all_outputs.append(channel_output)
  output =Image.merge('RGB', (all_outputs[0],all_outputs[1],all_outputs[2]))
  output.save(args[2])
  progress_bar.update(100-current_progress)
else:
  print("Invalid mode for input image:", mode)
  print("Only grayscale(L) and RGB(RGBA) images are supported.")


text:
"""
Note: includes tests for `last`
"""

import numpy as np
import pytest

from pandas import Series, date_range
import pandas._testing as tm

class TestFirst:
    def test_first_subset(self):
        rng = date_range("1/1/2000", "1/1/2010", freq="12h")
        ts = Series(np.random.randn(len(rng)), index=rng)
        result = ts.first("10d")
        assert len(result) == 20

        rng = date_range("1/1/2000", "1/1/2010", freq="D")
        ts = Series(np.random.randn(len(rng)), index=rng)
        result = ts.first("10d")
        assert len(result) == 10

        result = ts.first("3M")
        expected = ts[:"3/31/2000"]
        tm.assert_series_equal(result, expected)

        result = ts.first("21D")
        expected = ts[:21]
        tm.assert_series_equal(result, expected)

        result = ts[:0].first("3M")
        tm.assert_series_equal(result, ts[:0])

    def test_first_raises(self):
        # GH#20725
        ser = Series("a b c".split())
        msg = "'first' only supports a DatetimeIndex index"
        with pytest.raises(TypeError, match=msg):
            ser.first("1D")

    def test_last_subset(self):
        rng = date_range("1/1/2000", "1/1/2010", freq="12h")
        ts = Series(np.random.randn(len(rng)), index=rng)
        result = ts.last("10d")
        assert len(result) == 20

        rng = date_range("1/1/2000", "1/1/2010", freq="D")
        ts = Series(np.random.randn(len(rng)), index=rng)
        result = ts.last("10d")
        assert len(result) == 10

        result = ts.last("21D")
        expected = ts["12/12/2009":]
        tm.assert_series_equal(result, expected)

        result = ts.last("21D")
        expected = ts[-21:]
        tm.assert_series_equal(result, expected)

        result = ts[:0].last("3M")
        tm.assert_series_equal(result, ts[:0])

    def test_last_raises(self):
        # GH#20725
        ser = Series("a b c".split())
        msg = "'last' only supports a DatetimeIndex index"
        with pytest.raises(TypeError, match=msg):
            ser.last("1D")

max_stars_repo_licenses: MIT Unlicense

max_issues_repo_licenses: MIT Unlicense

max_forks_repo_licenses: MIT Unlicense