text:
import numpy as np

import medpy.metric.binary as mmb

import io_utils as io
import data_utils

import os

def evaluate_mmwhs(data_dir, data_list, model, label_ids, id_to_ignore = 0):
    # Evaluation function
    # Instead of re-training the model and reporting mean/std based on different random seeds, instead we just
    # report the mean/std based on the four test images. This seems to be the method used in
    # https://github.com/cchen-cc/SIFA/blob/master/evaluate.py
    # Thus, for each test image we obtain prediction labels, and find DICE and ASSM results
    # We make use of the pre-processing we did. Test tf-records 0..255 will correspond to image1,
    # test tf-records 256..511 will correspond to image2 etc.

    # Initialize result dicts
    dice = {}
    assd = {}
    for label in label_ids:
        if label_ids[label] == id_to_ignore:
            continue

        dice[label] = []
        assd[label] = []

    # Find results for each class
    for test_image in range(4):
        # The 3D labels and our predictions
        y_true = []
        y_hat = []

        for i in range(test_image * 256, (test_image+1) * 256):
            # We read in the i'th slice and compute predictions
            X_slice,y_true_slice = io.sample_batch(data_dir, [data_list[i]], 1)

            y_true.append(np.copy(y_true_slice))
            y_hat.append(model.predict(X_slice))

        y_true = np.array(y_true).reshape(256,256,256)
        y_hat = np.array(y_hat).reshape(256,256,256)

        print(y_true.shape, y_hat.shape)

        for label in label_ids:
            if label_ids[label] == id_to_ignore:
                continue

            # Prep data and compute metrics
            curr_y_true = np.copy(y_true)
            curr_y_true[curr_y_true != label_ids[label]] = 0

            curr_y_hat = np.copy(y_hat)
            curr_y_hat[curr_y_hat != label_ids[label]] = 0

            dice[label].append(mmb.dc(curr_y_hat, curr_y_true))
            assd[label].append(mmb.assd(curr_y_hat, curr_y_true))

    # Compute mean/std for each class for dice and assd
    dice_mu = {}
    dice_sd = {}
    dice_total_mu = []

    assd_mu = {}
    assd_sd = {}
    assd_total_mu = []

    for label in label_ids:
        if label_ids[label] == id_to_ignore:
            continue

        dice_mu[label] = np.mean(dice[label])
        dice_sd[label] = np.std(dice[label])
        dice_total_mu.append(dice_mu[label])

        assd_mu[label] = np.mean(assd[label])
        assd_sd[label] = np.std(assd[label])
        assd_total_mu.append(assd_mu[label])

    dice_total_mu = np.mean(dice_total_mu)
    assd_total_mu = np.mean(assd_total_mu)

    return dice_total_mu, dice_mu, dice_sd, assd_total_mu, assd_mu, assd_sd

def evaluate_abdomen(raw_ct_img_dir, data_seed, model, label_ids, id_to_ignore = 0, whole_dataset_computation=False):
    # Use the same random seed as in the data processing notebook
    assert data_seed == 0
    np.random.seed(data_seed)
    train_indices = sorted(np.random.choice(range(30), 24, replace=False))
    test_indices = np.asarray(sorted([x for x in range(30) if x not in train_indices]))
    test_images = np.asarray(sorted(os.listdir(raw_ct_img_dir + "img/")))[test_indices]
    test_labels = np.asarray(sorted(os.listdir(raw_ct_img_dir + "label/")))[test_indices]

    print(test_images)

    all_scans = []
    all_segmaps = []
    for imf,labelf in zip(test_images, test_labels):
        print("Working on a new scan", imf, labelf)

        scan,labels = data_utils.preprocess_abdomen_ct(raw_ct_img_dir, imf, labelf, final_cropping=True)

        all_scans.append(np.copy(scan))
        all_segmaps.append(np.copy(labels))

    # Compute results for all images at once
    if whole_dataset_computation == True:
        combined_scans = []
        for scan in all_scans:
            combined_scans.append(np.moveaxis(np.copy(scan), 2, 0))
        combined_scans = np.vstack(combined_scans)
        combined_scans = np.moveaxis(combined_scans, 0, 2)

        combined_segmaps = []
        for segmap in all_segmaps:
            combined_segmaps.append(np.moveaxis(np.copy(segmap), 2, 0))
        combined_segmaps = np.vstack(combined_segmaps)
        combined_segmaps = np.moveaxis(combined_segmaps, 0, 2)

        all_scans = [combined_scans]
        all_segmaps = [combined_segmaps]

    # Initialize result dicts
    dice = {}
    assd = {}
    for label in label_ids:
        if label_ids[label] == id_to_ignore:
            continue

        dice[label] = []
        assd[label] = []

    # Find results for each class
    for scan,labels in zip(all_scans, all_segmaps):
        # The 3D labels and our predictions
        y_true = []
        y_hat = []

        for idx in range(scan.shape[-1]):
            # Reshape the images/labels to 256x256xNum_Slices
            X,_ = io.get_consecutive_slices(scan, labels, idx, target_shape=(256,256,3))

            # Ignore slices that have none of the target organs present
            if len(np.unique(labels[...,idx])) == 1:
                continue

            Yhat = model.predict(X.reshape(1,256,256,3))[0]

            y_true.append(np.copy(labels[...,idx].reshape(1,256,256)))
            y_hat.append(np.copy(Yhat.reshape(1,256,256)))

        y_true = np.vstack(y_true)
        y_hat = np.vstack(y_hat)

        print(y_true.shape, y_hat.shape)

        for label in label_ids:
            if label_ids[label] == id_to_ignore:
                continue

            # Prep data and compute metrics
            curr_y_true = np.copy(y_true)
            curr_y_true[curr_y_true != label_ids[label]] = 0

            curr_y_hat = np.copy(y_hat)
            curr_y_hat[curr_y_hat != label_ids[label]] = 0

            dice[label].append(mmb.dc(curr_y_hat, curr_y_true))
            assd[label].append(mmb.assd(curr_y_hat, curr_y_true))

    # Compute mean/std for each class for dice and assd
    dice_mu = {}
    dice_sd = {}
    dice_total_mu = []

    assd_mu = {}
    assd_sd = {}
    assd_total_mu = []

    for label in label_ids:
        if label_ids[label] == id_to_ignore:
            continue

        dice_mu[label] = np.mean(dice[label])
        dice_sd[label] = np.std(dice[label])
        dice_total_mu.append(dice_mu[label])

        assd_mu[label] = np.mean(assd[label])
        assd_sd[label] = np.std(assd[label])
        assd_total_mu.append(assd_mu[label])

    dice_total_mu = np.mean(dice_total_mu)
    assd_total_mu = np.mean(assd_total_mu)

    return dice_total_mu, dice_mu, dice_sd, assd_total_mu, assd_mu, assd_sd


text:
import logging
import pickle
from collections import deque
from pathlib import Path

import h5py
import numpy as np

import Resources.training as r
from Utils.data_utils import extract_nearest_mesh_nodes_to_sensors, get_node_propery_at_states, \
    extract_coords_of_mesh_nodes
from Utils.img_utils import create_np_image

class DataloaderFlowfrontSensor:
    def __init__(self, image_size=None,
                 ignore_useless_states=True,
                 sensor_indizes=((0, 1), (0, 1)),
                 skip_indizes=(0, None, 1),
                 frame_count=1,
                 use_binary_sensor_only=False):
        self.image_size = image_size
        self.ignore_useless_states = ignore_useless_states
        self.sensor_indizes = sensor_indizes
        self.skip_indizes = skip_indizes
        self.mean = None
        self.std = None
        self.indeces_of_sensors = None
        self.frame_count = frame_count
        self.use_binary_sensor_only = use_binary_sensor_only
        self.coords = None
        self.load_ff_sensors()

    def load_ff_sensors(self):
        if r.nearest_nodes_to_sensors.is_file():
            with open(r.nearest_nodes_to_sensors, "rb") as nearest_nodes:
                _all_sensors = pickle.load(nearest_nodes)
        else:
            _all_sensors = extract_nearest_mesh_nodes_to_sensors(
                r.data_root / "2019-07-24_16-32-40_10p/0/2019-07-24_16-32-40_0")
            _all_sensors = _all_sensors.reshape((38, 30))
        with open(r.mean_std_20_flowfront_sensors, "rb") as mean_std:
            self.mean, self.std = pickle.load(mean_std)
        indices_of_sensors = _all_sensors[self.sensor_indizes[0][0]::self.sensor_indizes[0][1],
                                          self.sensor_indizes[1][0]::self.sensor_indizes[1][1]]
        self.indeces_of_sensors = indices_of_sensors.flatten()

    def extract_data_from_result_file(self, filename):
        with h5py.File(filename, 'r') as f:
            try:
                states = f["post"]["singlestate"]
                states = list(states)[self.skip_indizes[0]:self.skip_indizes[1]:self.skip_indizes[2]]

                filling_factors_at_certain_times = get_node_propery_at_states(f, "FILLING_FACTOR", states)
                velocity_at_certain_times = get_node_propery_at_states(f, "VELOCITY", states)
                fillings = np.squeeze(filling_factors_at_certain_times)
                velocities = np.squeeze(velocity_at_certain_times)
                self.extract_coords_data(f)
            except KeyError:
                logger = logging.getLogger()
                logger.warning(f'Warning: {filename}')
                return None, None, None
        return states, fillings, velocities

    def extract_data_from_meta_file(self, filename):
        with h5py.File(filename, 'r') as meta_file:
            try:
                useless_states = []
                if self.ignore_useless_states:
                    useless_states = meta_file["useless_states/singlestates"][()]
                array_of_states = meta_file["dryspot_states/singlestates"][()]
                set_of_dryspot_states = set(array_of_states.flatten())
            except KeyError:
                logger = logging.getLogger()
                logger.warning(f'Warning: {filename}')
                return None, None
        return useless_states, set_of_dryspot_states

    def extract_coords_data(self, f: h5py.File):
        if self.coords is not None:
            return self.coords
        self.coords = extract_coords_of_mesh_nodes(Path(f.filename))

    # TODO simplify code here, look at DataloaderImages for inspiration -> two generators (data, label) -> zip()
    def get_flowfront_sensor_and_flowfront_label(self, filename):
        """
        Load the flow front for the given states or all available states if states is None
        """
        states, fillings, velocities = self.extract_data_from_result_file(filename)
        meta_fn = str(filename).replace("RESULT.erfh5", "meta_data.hdf5")
        useless_states, set_of_dryspot_states = self.extract_data_from_meta_file(meta_fn)
        if states is None or \
                fillings is None or \
                velocities is None or \
                useless_states is None or \
                set_of_dryspot_states is None:
            return None
        instances = []
        frame_q = deque(maxlen=self.frame_count)
        for i, (velocity, filling, state) in enumerate(zip(velocities, fillings, states)):
            if self.ignore_useless_states \
                    and len(useless_states) > 0 \
                    and state == f'state{useless_states[0]:012d}':
                break
            binary_ff_sensor_values = np.ceil(filling[self.indeces_of_sensors])
            if self.use_binary_sensor_only:
                values = binary_ff_sensor_values
            else:
                values = self.add_velocity_factor(binary_ff_sensor_values, values, velocity)
            label = create_np_image(target_shape=self.image_size, norm_coords=self.coords, data=filling)
            if self.frame_count <= 1:
                instances.append((values, label))
            else:
                frame_q.append(values)
                """
                Stack the current frames of the queue so that each frame is in one channel and start
                using the data just after having enough data in the queue
                """
                if self.frame_count > i + 1:
                    continue
                else:
                    instances.append((np.stack(list(frame_q), axis=1), label))

        return instances

    def get_flowfront_sensor_bool_dryspot(self, filename):
        """
        Load the flow front for the given states or all available states if states is None
        """
        states, fillings, velocities = self.extract_data_from_result_file(filename)
        meta_fn = str(filename).replace("RESULT.erfh5", "meta_data.hdf5")
        useless_states, set_of_dryspot_states = self.extract_data_from_meta_file(meta_fn)
        if states is None or \
                fillings is None or \
                velocities is None or \
                useless_states is None or \
                set_of_dryspot_states is None:
            return None
        instances = []
        frame_q = deque(maxlen=self.frame_count)
        for i, (velocity, filling, state) in enumerate(zip(velocities, fillings, states)):
            if self.ignore_useless_states \
                    and len(useless_states) > 0 \
                    and state == f'state{useless_states[0]:012d}':
                break
            binary_ff_sensor_values = np.ceil(filling[self.indeces_of_sensors])
            if self.use_binary_sensor_only:
                values = binary_ff_sensor_values
            else:
                values = self.add_velocity_factor(binary_ff_sensor_values, velocity)
            label = 0
            if int(str(state).replace("state", "0")) in set_of_dryspot_states:
                label = 1
            if self.frame_count <= 1:
                instances.append((values, label))
            else:
                frame_q.append(values)
                """
                Stack the current frames of the queue so that each frame is in one channel and start
                using the data just after having enough data in the queue
                """
                if self.frame_count > i + 1:
                    continue
                else:
                    instances.append((np.stack(list(frame_q), axis=1), label))

        return instances

    def add_velocity_factor(self, binary_ff_sensor_values, velocity):
        velocities_at_sensors = velocity[self.indeces_of_sensors]
        abs_velocities_at_sensors = np.linalg.norm(velocities_at_sensors, axis=1)
        standardized_velocities = (abs_velocities_at_sensors - self.mean) / self.std
        return binary_ff_sensor_values * standardized_velocities

if __name__ == '__main__':
    dl = DataloaderFlowfrontSensor(sensor_indizes=((1, 8), (1, 8)))
    dl.get_flowfront_sensor_bool_dryspot(
        r'Y:\data\RTM\Leoben\sim_output\2019-07-24_16-32-40_10p\0\2019-07-24_16-32-40_0_RESULT.erfh5')

ext: py

max_stars_repo_name: isse-augsburg/rtm-predictions

max_issues_repo_name: isse-augsburg/rtm-predictions

max_forks_repo_name: isse-augsburg/rtm-predictions

reason: import numpy


text:
#!/usr/bin/env python
from __future__ import print_function
from __future__ import division
import base
import numpy as np

def tasks():
    return [
        AcidSubstitutions(),
    ]

class AcidSubstitutions(base.Task):
    """
    Creates data files for plots of the frequency with which amino acids are
    involved in substitutions, in our database.
    """
    def __init__(self):
        super(AcidSubstitutions, self).__init__('acid_substitutions')
        self._set_data_subdir('papergp')
        self._acids = None
        self._acid_map = None

    def acids(self):
        """
        Returns a list of amino acids and a reverse lookup dict.

        (So a list [A, R, N, ...] and a dict {A:0, R:1, N:2, ...}).
        """
        if self._acids is None:
            if False:
                acids = []
                acid_map = {}
                with base.connect() as con:
                    c = con.cursor()
                    q = 'select key from acid order by rowid'
                    for k, row in enumerate(c.execute(q)):
                        acids.append(row['key'])
                        acid_map[row['key']] = k
                self._acids = acids
                self._acid_map = acid_map
            else:
                self._acids = 'KDERHQNSPTGACWFLIMYV'
                assert(len(self._acids) == 20)
                self._acid_map = dict(
                    zip(self._acids, range(len(self._acids))))
        return self._acids, self._acid_map

    def expected(self, simple_weights=False):
        """
        Returns a 2d matrix of expected amino acid substitutions ratios.

        Uses weights based on single nucleotide substitutions.

        By default, these take into account rates derived from the human
        genome. To disable this, set simple_weights = True.
        """
        acids, acid_map = self.acids()
        n = len(acids)

        matrix = np.zeros((n, n))
        with base.connect() as con:
            c = con.cursor()
            q = 'select old, new, '
            q += 'simple_weight as weight' if simple_weights else 'weight'
            q += ' from mutation_possible'
            for k, row in enumerate(c.execute(q)):
                old = acid_map[row['old']]
                new = acid_map[row['new']]
                matrix[old, new] += row['weight']
                if old == new:
                    print(row['old'], row['new'])

        return matrix / np.sum(matrix)

    def observed(self):
        """
        Returns a matrix of observed amino acid substituion ratios.
        """
        acids, acid_map = self.acids()
        n = len(acids)

        matrix = np.zeros((n, n), dtype=float)
        with base.connect() as con:
            c = con.cursor()
            q = 'select distinct old, idx, new from report where pub != "exac"'
            for k, row in enumerate(c.execute(q)):
                old = acid_map[row['old']]
                new = acid_map[row['new']]
                matrix[old, new] += 1

        return matrix / np.sum(matrix)

    def calculate(self, basename, simple_weights):

        # Get list of acids, and a mapping from acid to list index
        acids, acid_map = self.acids()
        n = len(acids)

        # Get matrix of expected amino acid substitution ratios
        exp = self.expected()

        # Get matrix of observed amino acid substitution ratios
        obs = self.observed()

        # Store positions of
        #  - not observed, not expected
        #  - observed, not expected
        #  - expected, not observed
        not_obs_not_exp = (obs == 0) * (exp == 0)
        obs_not_exp = (obs != 0) * (exp == 0)
        exp_not_obs = (obs == 0) * (exp != 0)
        obs_exp = (obs != 0) * (exp != 0)
        assert(not np.any(not_obs_not_exp * obs_not_exp))
        assert(not np.any(not_obs_not_exp * exp_not_obs))
        assert(not np.any(obs_not_exp * exp_not_obs))

        # Calculate observed vs expected ratios
        with np.errstate(all='ignore'):
            ratios = obs / exp

        # Find lowest and highest (finite) ratio
        rmin = np.min(obs[obs_exp] / exp[obs_exp])
        rmax = np.max(obs[obs_exp] / exp[obs_exp])

        # Fill in special values for ratios
        ratios[not_obs_not_exp] = 1
        ratios[obs_not_exp] = rmax
        ratios[exp_not_obs] = rmin

        # Write file with ratios
        filename = self.data_out(basename + '.csv')
        print('Writing ' + filename)
        with open(filename, 'w') as f:
            w = self.csv_writer(f)
            for row in ratios:
                w.writerow(row)

        # Write exp-not-obs to file
        filename = self.data_out(basename + '-exp-not-obs.csv')
        print('Writing ' + filename)
        with open(filename, 'w') as f:
            w = self.csv_writer(f)
            w.writerow(('epx-not-obs-x', 'exp-not-obs-y'))
            for i in range(n):
                for j in range(n):
                    if exp_not_obs[i, j]:
                        w.writerow((0.5 + j, 19.5 - i))

        # Write obs-not-exp to file
        filename = self.data_out(basename + '-obs-not-exp.csv')
        print('Writing ' + filename)
        with open(filename, 'w') as f:
            w = self.csv_writer(f)
            w.writerow(('obs-not-exp-x', 'obs-not-exp-y'))
            for i in range(n):
                for j in range(n):
                    if obs_not_exp[i, j]:
                        w.writerow((0.5 + j, 19.5 - i))

        # Write not-obs-not-exp to file
        filename = self.data_out(basename + '-not-obs-not-exp.csv')
        print('Writing ' + filename)
        with open(filename, 'w') as f:
            w = self.csv_writer(f)
            w.writerow(('not-obs-not-exp-x', 'not-obs-not-exp-y'))

            for i in range(n):
                for j in range(n):
                    if not_obs_not_exp[i, j]:
                        w.writerow((0.5 + j, 19.5 - i))

        # Gather marginals
        # Note: No imputations needed here!
        # Indices in matrix are matrix[old, new]
        # So each row is for a single old acid, each column for a single new
        # Summing over axis=0 gives sum per column (new)
        # Summing over axis=1 gives sum per row (old)
        marginals_old = np.sum(obs, axis=1) / np.sum(exp, axis=1)
        marginals_new = np.sum(obs, axis=0) / np.sum(exp, axis=0)

        # Write file with margins
        filename = self.data_out(basename + '-marginals.csv')
        print('Writing ' + filename)
        # This also needs x and y-axis positions for the histograms
        xax = 0.5 + np.arange(20)
        yax = 19.5 - np.arange(20)
        with open(filename, 'w') as f:
            w = self.csv_writer(f)
            w.writerow((
                'marg-acid', 'marg-old', 'marg-new', 'marg-x', 'marg-y'))
            for i, acid in enumerate(acids):
                w.writerow((
                    acid, marginals_old[i], marginals_new[i], xax[i], yax[i]))

    def _run(self):

        self.calculate('substitutions-1-weighted', False)

if __name__ == '__main__':
    t = base.TaskRunner()
    t.add_tasks(tasks())
    t.run()


text:
# Copyright 2018 The Cirq Developers
#
# 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
#
#     https://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.

"""Abstract base classes for different types of simulators.

Simulator types include
    SimulatesSamples: mimics the interface of quantum hardware.
    SimulatesFinalWaveFunction: allows access to the wave function.
"""

import abc
import collections

from typing import Dict, Iterator, List, Union

import numpy as np

from cirq import circuits, ops, schedules, study
from cirq.sim import state

class SimulatesSamples:
    """Simulator that mimics running on quantum hardware.

    Implementors of this interface should implement the _run method.
    """

    def run(
        self,
        circuit: circuits.Circuit,
        param_resolver: study.ParamResolver = study.ParamResolver({}),
        repetitions: int = 1,
    ) -> study.TrialResult:
        """Runs the entire supplied Circuit, mimicking the quantum hardware.

        Args:
            circuit: The circuit to simulate.
            param_resolver: Parameters to run with the program.
            repetitions: The number of repetitions to simulate.

        Returns:
            TrialResult for a run.
        """
        return self.run_sweep(circuit, [param_resolver], repetitions)[0]

    def run_sweep(
        self,
        program: Union[circuits.Circuit, schedules.Schedule],
        params: study.Sweepable = study.ParamResolver({}),
        repetitions: int = 1,
    ) -> List[study.TrialResult]:
        """Runs the entire supplied Circuit, mimicking the quantum hardware.

        In contrast to run, this allows for sweeping over different parameter
        values.

        Args:
            program: The circuit or schedule to simulate.
            params: Parameters to run with the program.
            repetitions: The number of repetitions to simulate.

        Returns:
            TrialResult list for this run; one for each possible parameter
            resolver.
        """
        circuit = (program if isinstance(program, circuits.Circuit)
                   else program.to_circuit())
        param_resolvers = study.to_resolvers(params or study.ParamResolver({}))

        trial_results = []  # type: List[study.TrialResult]
        for param_resolver in param_resolvers:
            measurements = self._run(circuit=circuit,
                                     param_resolver=param_resolver,
                                     repetitions=repetitions)
            trial_results.append(study.TrialResult(params=param_resolver,
                                                   repetitions=repetitions,
                                                   measurements=measurements))
        return trial_results

    @abc.abstractmethod
    def _run(
        self,
        circuit: circuits.Circuit,
        param_resolver: study.ParamResolver,
        repetitions: int
    ) -> Dict[str, List[np.ndarray]]:
        """Run a simulation, mimicking quantum hardware.

        Args:
            circuit: The circuit to simulate.
            param_resolver: Parameters to run with the program.
            repetitions: Number of times to repeat the run.

        Returns:
            A dictionary from measurement key to a list of lists representing
            the results. Measurement results are a list of lists (a numpy
            ndarray), the first list corresponding to the repetition, and the
            second is the actual boolean measurement results (ordered by
            the qubits acted upon by the measurement gate.)
        """
        raise NotImplementedError()

class SimulatesFinalWaveFunction:
    """Simulator that allows access to a quantum computer's wavefunction.

    Implementors of this interface should implement the simulate_sweep
    method. This simulator only returns the wave function for the final
    step of a simulation. For simulators that also allow stepping through
    a circuit see `SimulatesIntermediateWaveFunction`.
    """

    def simulate(
        self,
        circuit: circuits.Circuit,
        param_resolver: study.ParamResolver = study.ParamResolver({}),
        qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT,
        initial_state: Union[int, np.ndarray] = 0,
    ) -> 'SimulationTrialResult':
        """Simulates the entire supplied Circuit.

        This method returns a result which allows access to the entire
        wave function.

        Args:
            circuit: The circuit to simulate.
            param_resolver: Parameters to run with the program.
            qubit_order: Determines the canonical ordering of the qubits used
                to define the order of amplitudes in the wave function.
            initial_state: If an int, the state is set to the computational
                basis state corresponding to this state. Otherwise  if this
                is a np.ndarray it is the full initial state. In this case it
                must be the correct size, be normalized (an L2 norm of 1), and
                be safely castable to an appropriate dtype for the simulator.

        Returns:
            SimulateTrialResults for the simulation. Includes the final wave
            function.
        """
        return self.simulate_sweep(circuit, [param_resolver], qubit_order,
                                   initial_state)[0]

    @abc.abstractmethod
    def simulate_sweep(
        self,
        program: Union[circuits.Circuit, schedules.Schedule],
        params: study.Sweepable = study.ParamResolver({}),
        qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT,
        initial_state: Union[int, np.ndarray] = 0,
    ) -> List['SimulationTrialResult']:
        """Simulates the entire supplied Circuit.

        This method returns a result which allows access to the entire
        wave function. In contrast to simulate, this allows for sweeping
        over different parameter values.

        Args:
            program: The circuit or schedule to simulate.
            params: Parameters to run with the program.
            qubit_order: Determines the canonical ordering of the qubits used to
                define the order of amplitudes in the wave function.
            initial_state: If an int, the state is set to the computational
                basis state corresponding to this state.
                Otherwise if this is a np.ndarray it is the full initial state.
                In this case it must be the correct size, be normalized (an L2
                norm of 1), and  be safely castable to an appropriate
                dtype for the simulator.

        Returns:
            List of SimulatorTrialResults for this run, one for each
            possible parameter resolver.
        """
        raise NotImplementedError()

class SimulationTrialResult:
    """Results of a simulation by a SimulatesFinalWaveFunction.

    Unlike TrialResult these results contain the final state (wave function)
    of the system.

    Attributes:
        params: A ParamResolver of settings used for this result.
        measurements: A dictionary from measurement gate key to measurement
            results. Measurement results are a numpy ndarray of actual boolean
            measurement results (ordered by the qubits acted on by the
            measurement gate.)
        final_state: The final state (wave function) of the system after the
            trial finishes.
    """

    def __init__(self,
        params: study.ParamResolver,
        measurements: Dict[str, np.ndarray],
        final_state: np.ndarray) -> None:
        self.params = params
        self.measurements = measurements
        self.final_state = final_state

    def __repr__(self):
        return ('SimulationTrialResult(params={!r}, '
                'measurements={!r}, '
                'final_state={!r})').format(self.params,
                                            self.measurements,
                                            self.final_state)

    def __str__(self):
        def bitstring(vals):
            return ''.join('1' if v else '0' for v in vals)

        results = sorted(
            [(key, bitstring(val)) for key, val in self.measurements.items()])
        return ' '.join(
            ['{}={}'.format(key, val) for key, val in results])

    def dirac_notation(self, decimals=2):
        return state.dirac_notation(self.final_state, decimals)

    def _eq_tuple(self):
        measurements = {k: v.tolist() for k, v in
                        sorted(self.measurements.items())}
        return (SimulationTrialResult, self.params, measurements,
                self.final_state.tolist())

    def __eq__(self, other):
        if not isinstance(other, SimulationTrialResult):
            return NotImplemented
        return self._eq_tuple() == other._eq_tuple()

    def __ne__(self, other):
        return not self == other

class SimulatesIntermediateWaveFunction(SimulatesFinalWaveFunction):
    """A SimulatesFinalWaveFunction that simulates a circuit by moments.

    Whereas a general SimulatesFinalWaveFunction may return the entire wave
    function at the end of a circuit, a SimulatesIntermediateWaveFunction can
    simulate stepping through the moments of a circuit.

    Implementors of this interface should implement the _simulator_iterator
    method.
    """

    def simulate_sweep(
        self,
        program: Union[circuits.Circuit, schedules.Schedule],
        params: study.Sweepable = study.ParamResolver({}),
        qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT,
        initial_state: Union[int, np.ndarray] = 0,
    ) -> List['SimulationTrialResult']:
        """Simulates the entire supplied Circuit.

        This method returns a result which allows access to the entire
        wave function. In contrast to simulate, this allows for sweeping
        over different parameter values.

        Args:
            program: The circuit or schedule to simulate.
            params: Parameters to run with the program.
            qubit_order: Determines the canonical ordering of the qubits used to
                define the order of amplitudes in the wave function.
            initial_state: If an int, the state is set to the computational
                basis state corresponding to this state.
                Otherwise if this is a np.ndarray it is the full initial state.
                In this case it must be the correct size, be normalized (an L2
                norm of 1), and  be safely castable to an appropriate
                dtype for the simulator.

        Returns:
            List of SimulatorTrialResults for this run, one for each
            possible parameter resolver.
        """
        circuit = (program if isinstance(program, circuits.Circuit)
                   else program.to_circuit())
        param_resolvers = study.to_resolvers(params or study.ParamResolver({}))

        trial_results = []  # type: List[SimulationTrialResult]
        qubit_order = ops.QubitOrder.as_qubit_order(qubit_order)
        for param_resolver in param_resolvers:
            step_result = None
            all_step_results = self.simulate_moment_steps(circuit,
                                                          param_resolver,
                                                          qubit_order,
                                                          initial_state)
            measurements = {}  # type: Dict[str, np.ndarray]
            for step_result in all_step_results:
                for k, v in step_result.measurements.items():
                    measurements[k] = np.array(v, dtype=bool)
            if step_result:
                final_state = step_result.state()
            else:
                # Empty circuit, so final state should be initial state.
                print(circuit)
                num_qubits = len(qubit_order.order_for(circuit.all_qubits()))
                final_state = state.to_valid_state_vector(initial_state,
                                                          num_qubits)
            trial_results.append(SimulationTrialResult(
                params=param_resolver,
                measurements=measurements,
                final_state=final_state))

        return trial_results

    def simulate_moment_steps(
        self,
        circuit: circuits.Circuit,
        param_resolver: study.ParamResolver = None,
        qubit_order: ops.QubitOrderOrList = ops.QubitOrder.DEFAULT,
        initial_state: Union[int, np.ndarray] = 0
    ) -> Iterator['StepResult']:
        """Returns an iterator of StepResults for each moment simulated.

        Args:
            circuit: The Circuit to simulate.
            param_resolver: A ParamResolver for determining values of Symbols.
            qubit_order: Determines the canonical ordering of the qubits used to
                define the order of amplitudes in the wave function.
            initial_state: If an int, the state is set to the computational
                basis state corresponding to this state. Otherwise if this is
                a np.ndarray it is the full initial state. In this case it must
                be the correct size, be normalized (an L2 norm of 1), and
                be safely castable to an appropriate dtype for the simulator.

        Returns:
            Iterator that steps through the simulation, simulating each
            moment and returning a StepResult for each moment.
        """
        param_resolver = param_resolver or study.ParamResolver({})
        return self._simulator_iterator(circuit, param_resolver, qubit_order,
                                        initial_state)

    @abc.abstractmethod
    def _simulator_iterator(
        self,
        circuit: circuits.Circuit,
        param_resolver: study.ParamResolver,
        qubit_order: ops.QubitOrderOrList,
        initial_state: Union[int, np.ndarray],
    ) -> Iterator['StepResult']:
        """Iterator over StepResult from Moments of a Circuit.

        Args:
            circuit: The circuit to simulate.
            param_resolver: A ParamResolver for determining values of
                Symbols.
            qubit_order: Determines the canonical ordering of the qubits used to
                define the order of amplitudes in the wave function.
            initial_state: The full initial state. This must be the correct
                size, be normalized (an L2 norm of 1), and be safely
                castable to a complex type handled by the simulator.

        Yields:
            StepResults from simulating a Moment of the Circuit.
        """
        raise NotImplementedError()

class StepResult:
    """Results of a step of a SimulatesFinalWaveFunction.

    Attributes:
        qubit_map: A map from the Qubits in the Circuit to the the index
            of this qubit for a canonical ordering. This canonical ordering is
            used to define the state (see the state() method).
        measurements: A dictionary from measurement gate key to measurement
            results, ordered by the qubits that the measurement operates on.
    """

    def __init__(
        self,
        qubit_map: Dict,
        measurements: Dict[str, List[bool]]) -> None:
        self.qubit_map = qubit_map or {}
        self.measurements = measurements or collections.defaultdict(list)

    @abc.abstractmethod
    def state(self) -> np.ndarray:
        """Return the state (wave function) at this point in the computation.

        The state is returned in the computational basis with these basis
        states defined by the qubit_map. In particular the value in the
        qubit_map is the index of the qubit, and these are translated into
        binary vectors where the last qubit is the 1s bit of the index, the
        second-to-last is the 2s bit of the index, and so forth (i.e. big
        endian ordering).

        Example:
             qubit_map: {QubitA: 0, QubitB: 1, QubitC: 2}
             Then the returned vector will have indices mapped to qubit basis
             states like the following table
               |   | QubitA | QubitB | QubitC |
               +---+--------+--------+--------+
               | 0 |   0    |   0    |   0    |
               | 1 |   0    |   0    |   1    |
               | 2 |   0    |   1    |   0    |
               | 3 |   0    |   1    |   1    |
               | 4 |   1    |   0    |   0    |
               | 5 |   1    |   0    |   1    |
               | 6 |   1    |   1    |   0    |
               | 7 |   1    |   1    |   1    |
               +---+--------+--------+--------+
        """
        raise NotImplementedError()

    @abc.abstractmethod
    def set_state(self, state: Union[int, np.ndarray]):
        """Updates the state of the simulator to the given new state.

        Args:
            state: If this is an int, then this is the state to reset
            the stepper to, expressed as an integer of the computational basis.
            Integer to bitwise indices is little endian. Otherwise if this is
            a np.ndarray this must be the correct size and have dtype of
            np.complex64.

        Raises:
            ValueError if the state is incorrectly sized or not of the correct
            dtype.
        """
        raise NotImplementedError()

    @abc.abstractmethod
    def sample(self, qubits: List[ops.QubitId], repetitions: int = 1):
        """Samples from the wave function at this point in the computation.

        Note that this does not collapse the wave function.

        Returns:
            Measurement results with True corresponding to the |1> state.
            The outer list is for repetitions, and the inner corresponds to
            measurements ordered by the supplied qubits.
        """
        raise NotImplementedError()

    def dirac_notation(self, decimals=2):
        """Returns the wavefunction as a string in Dirac notation.

        Args:
            state: A sequence representing a wave function in which the ordering
                mapping to qubits follows the standard Kronecker convention of
                numpy.kron.
            decimals: How many decimals to include in the pretty print.

        Returns:
            A pretty string consisting of a sum of computational basis kets
            and non-zero floats of the specified accuracy."""
        return state.dirac_notation(self.state, decimals)

max_stars_repo_name: sboixo/Cirq

max_issues_repo_name: mannurulz/Cirq

max_forks_repo_name: mannurulz/Cirq


text:
from __future__ import print_function, division

import os
import shutil
import tempfile

import pytest
import numpy as np
from numpy.testing import assert_array_almost_equal_nulp

from astropy import units as u

from .. import Model
from ..image import Image
from ...util.functions import random_id
from .test_helpers import get_test_dust, get_highly_reflective_dust

class TestFilters(object):

    def setup_class(self):

        m = Model()

        m.set_cartesian_grid([-1., 1.],
                             [-1., 1.],
                             [-1., 1.])

        m.add_density_grid(np.array([[[1.]]]), get_test_dust())

        s = m.add_point_source()
        s.name = 'first'
        s.luminosity = 1.
        s.temperature = 6000.

        s = m.add_point_source()
        s.name = 'second'
        s.luminosity = 1.
        s.temperature = 6000.

        i = m.add_peeled_images(sed=True, image=True)
        i.set_viewing_angles([1., 2., 3.], [1., 2., 3.])
        i.set_image_limits(-1., 1., -1., 1.)
        i.set_image_size(10, 20)

        f1 = i.add_filter()
        f1.name = 'F1'
        f1.spectral_coord = [1, 1.1, 1.2, 1.3] * u.micron
        f1.transmission = [0., 100., 50, 0.] * u.percent
        f1.detector_type = 'photons'
        f1.alpha = 0.
        f1.central_spectral_coord = 1.15 * u.micron

        f2 = i.add_filter()
        f2.name = 'F2'
        f2.spectral_coord = [2, 2.1, 2.2, 2.3, 2.4] * u.micron
        f2.transmission = [0., 50, 100, 60, 0.] * u.percent
        f2.detector_type = 'energy'
        f2.alpha = 1.
        f2.central_spectral_coord = 2.15 * u.micron

        m.set_n_initial_iterations(0)

        m.set_n_photons(imaging=1000)

        self.tmpdir = tempfile.mkdtemp()
        m.write(os.path.join(self.tmpdir, random_id()))

        self.m = m.run()

    def teardown_class(self):
        shutil.rmtree(self.tmpdir)

    def test_image_wav(self):
        image = self.m.get_image()
        np.testing.assert_allclose(image.nu, [2.60689094e+14, 1.39438353e+14])
        np.testing.assert_allclose(image.wav, [1.15, 2.15])

    def test_sed_wav(self):
        sed = self.m.get_sed()
        np.testing.assert_allclose(sed.nu, [2.60689094e+14, 1.39438353e+14])
        np.testing.assert_allclose(sed.wav, [1.15, 2.15])

    def test_image_shape(self):
        image = self.m.get_image()
        assert image.val.shape == (3, 20, 10, 2)

    def test_sed_shape(self):
        sed = self.m.get_sed()
        assert sed.val.shape == (3, 1, 2)

    def test_image_values(self):
        image = self.m.get_image(units='MJy/sr', distance=1)
        np.testing.assert_allclose(np.sum(image.val[:, :, :, 0]), 3438.059082285024, rtol=0.1)
        np.testing.assert_allclose(np.sum(image.val[:, :, :, 1]), 2396.4803378036186, rtol=0.1)

max_stars_repo_name: christopherlovell/hyperion

max_issues_repo_name: christopherlovell/hyperion

max_forks_repo_name: christopherlovell/hyperion

reason: import numpy,from numpy,from astropy


text:
import os
import numpy as np
import tensorflow as tf
# tf.logging.set_verbosity(tf.logging.INFO)
import pdb
from tqdm import trange
from parse import load_data

class FRRN():

    def __init__(self, x_shp, model='A'):

        self.Training = True
        self.x_shp = x_shp
        self.num_class = 12 # we will have an additional label for negative
        self.reg_lambda = 1e-4
        self.learning_rate = 1e-4
        self.train_batch_size = 1
        # if valid batch size too big, GPU boom
        self.valid_batch_size = 1
        self.test_batch_size = 1
        self.validation_frequency = 10

        self.report_freq = 5
        # top K pixels loss
        self.K = 512 * 16
        self.max_iter = 50

        self.log_dir = 'logs'
        self.save_dir = 'save'

        # build the network
        self._build_placeholder()
        self._build_preprocessing()
        self._build_model_a() if model == 'A' else self._build_model_b()
        self._build_loss()
        self._build_optim()
        self._build_eval()
        self._build_summary()
        self._build_writer()

    def _build_placeholder(self):
        # Build placeholders.

        # N * H * W * C
        x_in_shp = (self.train_batch_size, self.x_shp[1], self.x_shp[2], self.x_shp[3])

        # Create Placeholders for inputs
        self.x_in = tf.placeholder(tf.float32, shape=x_in_shp)
        self.y_in = tf.placeholder(tf.int64, shape=x_in_shp[:3])

    def _build_preprocessing(self):
        """Build preprocessing related graph."""

        with tf.variable_scope("Normalization", reuse=tf.AUTO_REUSE):
            # TODO: we will make `n_mean`, `n_range`, `n_mean_in` and
            # `n_range_in` as scalar this time! This is how we often use in
            # CNNs, as we KNOW that these are image pixels, and all pixels
            # should be treated equally!

            # Create placeholders for saving mean, range to a TF variable for
            # easy save/load. Create these variables as well.
            self.n_mean_in = tf.placeholder(tf.float32, shape=())
            self.n_range_in = tf.placeholder(tf.float32, shape=())
            # Make the normalization as a TensorFlow variable. This is to make
            # sure we save it in the graph
            self.n_mean = tf.get_variable(
                "n_mean", shape=(), trainable=False)
            self.n_range = tf.get_variable(
                "n_range", shape=(), trainable=False)
            # Assign op to store this value to TF variable
            self.n_assign_op = tf.group(
                tf.assign(self.n_mean, self.n_mean_in),
                tf.assign(self.n_range, self.n_range_in),
            )

    # --z------residual stream->---------z_prime--->
    # 	|							^
    # 	|							|
    # 	V							|
    # --y------pooling stream-->---------y_prime--->

    def frr_unit(self, residual_stream, pooling_stream, max_pooling_pool_size, max_pooling_strides,
                max_pooling_padding, conv2d_num_filter, conv2d_kernel, conv2d_strides, conv2d_padding, activ,
                training=True):
        # iter1: pooling_stream[?, 360, 480, 48], residual stream[?, 720, 960, 32]-->[?, 360, 480, 32] by maxpooling
        residual_stream_orig = residual_stream
        residual_stream = tf.layers.max_pooling2d(residual_stream, max_pooling_pool_size, max_pooling_strides,
                         padding=max_pooling_padding)
        #concatenate 2 stream for conv2d
        # iter1: pooling_stream[?, 360, 480, 80]
        pooling_stream = tf.concat([pooling_stream, residual_stream], axis = 3)
        # iter1: pooling_stream[?, 358, 478, 96]
        pooling_stream = activ(tf.layers.conv2d(pooling_stream, conv2d_num_filter, conv2d_kernel, strides=conv2d_strides,
                         padding=conv2d_padding))
        pooling_stream = tf.layers.batch_normalization(pooling_stream, axis=-1, momentum=1.0, center=True, scale=True,
                         training=training, trainable=True, reuse=None)
        # iter1: pooling_stream[?, 356, 476, 96]
        pooling_stream = activ(tf.layers.conv2d(pooling_stream, conv2d_num_filter, conv2d_kernel, strides=conv2d_strides,
                         padding=conv2d_padding))
        pooling_stream = tf.layers.batch_normalization(pooling_stream, axis=-1, momentum=1.0, center=True, scale=True,
                         training=training, trainable=True, reuse=None)
        y_prime = tf.nn.relu(pooling_stream)
        # this one have same kernel and strids
        # iter1: pooling_stream[?, 356, 476, 32]
        pooling_stream = activ(tf.layers.conv2d(y_prime, 32, (1, 1), strides=(1, 1), padding='same'))
        # iter1: pooling_stream[?, 712, 952, 32]
        pooling_stream = tf.image.resize_nearest_neighbor(pooling_stream, tf.shape(y_prime)[1:3] * max_pooling_strides)
        # concatenate
        z_prime = residual_stream_orig + pooling_stream

        return z_prime, y_prime

    def residual_unit(self, fan_in, num_filter, kernel, strides, padding, activ, traning=True):

        cur_in = activ(tf.layers.conv2d(fan_in, num_filter, kernel, strides=strides, padding=padding))
        cur_in = tf.layers.batch_normalization(cur_in, axis=-1, momentum=1.0, center=True, scale=True, training=traning,
                                                 trainable=True, reuse=None)
        cur_in = tf.nn.relu(cur_in)
        cur_in = activ(tf.layers.conv2d(cur_in, num_filter, kernel, strides=strides, padding=padding))
        cur_in = tf.layers.batch_normalization(cur_in, axis=-1, momentum=1.0, center=True, scale=True, training=traning,
                                                 trainable=True, reuse=None)

        return cur_in + fan_in

    def _build_model_a(self):

        activ = tf.nn.relu
        kernel_initializer = tf.keras.initializers.he_normal()

        with tf.variable_scope("FRRN_A", reuse = tf.AUTO_REUSE):
            # batch_norm
            cur_in = (self.x_in - self.n_mean) / self.n_range
            # before 3 * ResNet
            cur_in = tf.layers.conv2d(cur_in, 48, (5, 5), strides=(1, 1), padding='same', name='conv2d_1')
            cur_in = tf.layers.batch_normalization(cur_in, axis=-1, momentum=1.0, center=True,
                     scale=True, training=self.Training, trainable=True, name="conv2d_1_bn", reuse=None)
            cur_in = activ(cur_in)
            # 3 * RU
            for i in range(3):
                cur_in = self.residual_unit(cur_in, 48, (3, 3), (1, 1), 'same', activ, traning=True)

            # now the shape should be as [?, 720, 960, 48]->[?, 720, 960, 32]. Splited to cur_in & residual_stream
            residual_stream = activ(tf.layers.conv2d(cur_in, 32, (1, 1), strides=(1, 1), padding='same',
                     name='conv2d_to_split'))

            # encoding
            pooling_stream = cur_in
            for it, num_filter, scale in [(3, 96, 2), (4, 192, 4), (2, 384, 8), (2, 384, 16)]:
                # max pool the pooling stream only
                # iter1: [?, 720, 960, 48]->[?, 360, 480, 48]
                pooling_stream = tf.layers.max_pooling2d(pooling_stream, (2, 2), (2, 2), padding='same')

                for i in range(it):
                    residual_stream, pooling_stream = self.frr_unit(residual_stream, pooling_stream, scale, scale, 'same',
                                                            num_filter, (3, 3), (1, 1), 'same', activ, training=self.Training)

            # decoding
            for it, num_filter, scale in [(2, 192, 8),(2, 192, 4),(2, 96, 2)]:
                # invert the max pool operation
                pooling_stream = tf.image.resize_nearest_neighbor(pooling_stream, tf.shape(pooling_stream)[1:3] * 2)
                # in the very last iteration, we discard pooing_stream actually, back to residual only
                for i in range(it):
                    residual_stream, pooling_stream = self.frr_unit(residual_stream, pooling_stream, scale, scale, 'same',
                                                            num_filter, (3, 3), (1, 1), 'same', activ, training=self.Training)

            # encoding pooled 4 times, 3 times in decoding, so do the last unpooling
            residual_stream = tf.image.resize_nearest_neighbor(residual_stream, tf.shape(pooling_stream)[1:3] * 2)
            # Concat Streams, make the residual stream back to 48 channels, just rename it back to cur_in
            cur_in = activ(tf.layers.conv2d(residual_stream, 48, (1, 1), strides=(1, 1), padding='same', name='conv2d_to_merge'))
            # 3 * RU again
            for i in range(3):
                cur_in = self.residual_unit(cur_in, 48, (3, 3), (1, 1), 'same', activ, traning=True)

            # Final classification layer
            self.logits = activ(tf.layers.conv2d(cur_in, self.num_class + 1, (1, 1), strides=(1, 1), padding='same',
                             name='output_layer'))

            # Get list of all weights in this scope. They are called "kernel" in tf.layers.dense.
            self.kernels_list = [_v for _v in tf.trainable_variables() if "kernel" in _v.name]

    def _build_model_b(self):

        activ = tf.nn.relu
        kernel_initializer = tf.keras.initializers.he_normal()

        with tf.variable_scope("FRRN_B", reuse = tf.AUTO_REUSE):
            # batch_norm
            cur_in = (self.x_in - self.n_mean) / self.n_range
            # before 3 * ResNet
            cur_in = tf.layers.conv2d(cur_in, 48, (5, 5), strides=(1, 1), padding='same', name='conv2d_1')
            cur_in = tf.layers.batch_normalization(cur_in, axis=-1, momentum=1.0, center=True, scale=True, training=self.Training,
                                                     trainable=True, name="conv2d_1_bn", reuse=None)
            cur_in = activ(cur_in)
            # 3 * RU
            for i in range(3):
                cur_in = self.residual_unit(cur_in, 48, (3, 3), (1, 1), 'same', activ, traning=True)

            # now the shape should be as [?, 720, 960, 48]->[?, 720, 960, 32]. Splited to cur_in & residual_stream
            residual_stream = activ(tf.layers.conv2d(cur_in, 32, (1, 1), strides=(1, 1), padding='same', name='conv2d_to_split'))

            # encoding
            pooling_stream = cur_in
            for it, num_filter, scale in [(3, 96, 2), (4, 192, 4), (2, 384, 8), (2, 384, 16), (2, 384, 32)]:
                # max pool the pooling stream only
                # iter1: [?, 720, 960, 48]->[?, 360, 480, 48]
                pooling_stream = tf.layers.max_pooling2d(pooling_stream, (2, 2), (2, 2), padding='same')

                for i in range(it):
                    residual_stream, pooling_stream = self.frr_unit(residual_stream, pooling_stream, scale, scale, 'same',
                                                                     num_filter, (3, 3), (1, 1), 'same', activ, training=self.Training)

            # decoding
            for it, num_filter, scale in [(2, 192, 16),(2, 192, 8),(2, 192, 4),(2, 96, 2)]:
                # invert the max pool operation
                pooling_stream = tf.image.resize_nearest_neighbor(pooling_stream, tf.shape(pooling_stream)[1:3] * 2)
                # in the very last iteration, we discard pooing_stream actually, back to residual only
                for i in range(it):
                    residual_stream, pooling_stream = self.frr_unit(residual_stream, pooling_stream, scale, scale, 'same',
                                                                     num_filter, (3, 3), (1, 1), 'same', activ, training=self.Training)

            # encoding pooled 4 times, 3 times in decoding, so do the last unpooling
            residual_stream = tf.image.resize_nearest_neighbor(residual_stream, tf.shape(pooling_stream)[1:3] * 2)
            # Concat Streams, make the residual stream back to 48 channels, just rename it back to cur_in
            cur_in = activ(tf.layers.conv2d(residual_stream, 48, (1, 1), strides=(1, 1), padding='same', name='conv2d_to_merge'))
            # 3 * RU again
            for i in range(3):
                cur_in = self.residual_unit(cur_in, 48, (3, 3), (1, 1), 'same', activ, traning=True)

            # Final classification layer
            self.logits = activ(tf.layers.conv2d(cur_in, self.num_class + 1, (1, 1), strides=(1, 1), padding='same', name='output_layer'))

            # Get list of all weights in this scope. They are called "kernel" in tf.layers.dense.
            self.kernels_list = [_v for _v in tf.trainable_variables() if "kernel" in _v.name]

    def _build_loss(self):
        """Build our cross entropy loss."""

        with tf.variable_scope("Loss", reuse=tf.AUTO_REUSE):
            # Create cross entropy loss, shape of naive_loss is N * H * W * 1
            self.logits = tf.check_numerics(self.logits, "logits contains Nan !!!!!")
            self.naive_loss = tf.nn.sparse_softmax_cross_entropy_with_logits(labels=self.y_in, logits=self.logits)
            self.naive_loss = tf.Print(self.naive_loss, [tf.reduce_min(self.naive_loss)], "naive_loss_min")
            self.naive_loss = tf.Print(self.naive_loss, [tf.reduce_max(self.naive_loss)], "naive_loss_max")
            self.naive_loss = tf.check_numerics(self.naive_loss, "naive_loss contains Nan !!!!!")
            self.loss_shape = (tf.shape(self.logits)[0], tf.shape(self.logits)[1] * tf.shape(self.logits)[2])
            # boot strapped corss-entropy loss
            self.loss_reshaped, _ = tf.nn.top_k(tf.reshape(self.naive_loss, self.loss_shape), self.K, sorted=False)
            self.loss = tf.reduce_mean(self.naive_loss)
            # self.loss = tf.Print(self.loss, [self.loss], "loss")

            # Create l2 regularizer loss and add
            l2_loss = tf.add_n([tf.reduce_sum(_v**2) for _v in self.kernels_list])
            self.loss += self.reg_lambda * l2_loss

            # Record summary for loss
            tf.summary.scalar("loss", self.loss)

    def _build_optim(self):
        """Build optimizer related ops and vars."""

        with tf.variable_scope("Optim", reuse=tf.AUTO_REUSE):
            self.global_step = tf.get_variable("global_step", shape=(), initializer=tf.zeros_initializer(), dtype=tf.int64, trainable=False)
            self.optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            self.optim = self.optimizer.minimize(self.loss, global_step=self.global_step)

    def _build_eval(self):
        """Build the evaluation related ops"""

        with tf.variable_scope("Eval", tf.AUTO_REUSE):

            # Compute the accuracy of the model. When comparing labels
            # elemwise, use tf.equal instead of `==`. `==` will evaluate if
            # your Ops are identical Ops.
            self.pred = tf.argmax(self.logits, axis=3)
            self.acc = tf.reduce_mean(
                tf.to_float(tf.equal(self.pred, self.y_in))
            )

            # Record summary for accuracy
            tf.summary.scalar("accuracy", self.acc)

            # We also want to save best validation accuracy. So we do
            # something similar to what we did before with n_mean. Note that
            # these will also be a scalar variable
            self.best_va_acc_in = tf.placeholder(tf.float32, shape=())
            self.best_va_acc = tf.get_variable("best_va_acc", shape=(), trainable=False)
            self.acc_assign_op = tf.assign(self.best_va_acc, self.best_va_acc_in)

    def _build_summary(self):
        """Build summary ops."""

        # Merge all summary op
        self.summary_op = tf.summary.merge_all()

    def _build_writer(self):
        """Build the writers and savers"""

        # Create summary writers (one for train, one for validation)
        self.summary_tr = tf.summary.FileWriter(os.path.join(self.log_dir, "train"))
        self.summary_va = tf.summary.FileWriter(os.path.join(self.log_dir, "valid"))
        # Create savers (one for current, one for best)
        self.saver_cur = tf.train.Saver()
        self.saver_best = tf.train.Saver()
        # Save file for the current model
        self.save_file_cur = os.path.join(
            self.log_dir, "model")
        # Save file for the best model
        self.save_file_best = os.path.join(
            self.save_dir, "model")

    def train(self, x_tr, y_tr, x_va, y_va):

        x_tr_mean = x_tr.mean()
        x_tr_range = 128.0

        # ----------------------------------------
        # Run TensorFlow Session
        with tf.Session() as sess:
            # Init
            print("Initializing...")
            sess.run(tf.global_variables_initializer())

            # Assign normalization variables from statistics of the train data
            sess.run(self.n_assign_op, feed_dict={
                self.n_mean_in: x_tr_mean,
                self.n_range_in: x_tr_range,
            })

            # Check if previous train exists
            b_resume = tf.train.latest_checkpoint(self.log_dir)
            if b_resume:
                # Restore network
                print("Restoring from {}...".format(
                    self.log_dir))
                self.saver_cur.restore(
                    sess,
                    b_resume
                    )
                # restore number of steps so far
                step = sess.run(self.global_step)
                # restore best acc
                best_acc = sess.run(self.best_va_acc)
            else:
                print("Starting from scratch...")
                step = 0
                best_acc = 0

            print("Training...")
            # for each epoch
            for step in trange(step, self.max_iter):

                self.training=True
                # Get a random training batch. Notice that we are now going to
                # forget about the `epoch` thing. Theoretically, they should do
                # almost the same.
                ind_cur = np.random.choice(
                    len(x_tr), self.train_batch_size, replace=False)
                x_b = np.array([x_tr[_i] for _i in ind_cur])
                y_b = np.array([y_tr[_i] for _i in ind_cur])
                y_b[y_b < 0] = self.num_class

                b_write_summary = ((step + 1) % self.report_freq == 0) or (step == 0)

                if b_write_summary:
                    fetches = {
                        "optim": self.optim,
                        "loss": self.loss,
                        "summary": self.summary_op,
                        "global_step": self.global_step,
                    }
                else:
                    fetches = {
                        "optim": self.optim,
                        "loss": self.loss,
                    }

                # Run the operations necessary for training
                res = sess.run(
                    fetches=fetches,
                    feed_dict={
                        self.x_in: x_b,
                        self.y_in: y_b,
                    },
                )
                print("iteration {} training loss: {}".format(step, res['loss']))
                #
                if 'global_step' in res:
                    print('global_step is: ', res['global_step'])

                    self.summary_tr.add_summary(res["summary"], global_step=res["global_step"],)
                    self.summary_tr.flush()

                    # Also save current model to resume when we write the
                    # summary.
                    self.saver_cur.save(
                        sess, self.save_file_cur,
                        global_step=self.global_step,
                        write_meta_graph=False,
                    )

                b_validate = ((step + 1) % self.validation_frequency == 0) or (step == 0)
                if b_validate:

                    self.traning = False
                    num_valid_b = len(x_va) // self.valid_batch_size
                    acc_list = []
                    loss_list = []

                    for idx_b in range(num_valid_b):
                        res = sess.run(
                            fetches={
                                "acc": self.acc,
                                "loss": self.loss,
                                "global_step": self.global_step,
                            },
                            feed_dict={
                                self.x_in: x_va[idx_b * self.valid_batch_size: (idx_b + 1) * self.valid_batch_size],
                                self.y_in: y_va[idx_b * self.valid_batch_size: (idx_b + 1) * self.valid_batch_size],
                            }
                        )
                        print('valid batch avg loss: ', res['loss'])
                        acc_list += [res['acc']]
                        loss_list += [res['loss']]
                        # pdb.set_trace()
                    # Write Validation Summary over batches
                    avg_loss = np.mean(loss_list)
                    avg_acc = np.mean(acc_list)
                    summary = tf.Summary(value=[
                    	tf.Summary.Value(tag="accuracy", simple_value=avg_acc),
                    	tf.Summary.Value(tag="loss", simple_value=avg_loss),
                    	])
                    self.summary_va.add_summary(summary, global_step=res["global_step"],)
                    self.summary_va.flush()

                    # If best validation accuracy, update W_best, b_best, and
                    # best accuracy. We will only return the best W and b

                    if avg_acc > best_acc:
                        best_acc = avg_acc
                        sess.run(self.acc_assign_op,
                            feed_dict = {
                                self.best_va_acc_in: best_acc,
                            })
                        # Save the best model
                        self.saver_best.save(
                            sess, self.save_file_best,
                            write_meta_graph=False,
                            )

            print("best training validation accuracy is: {}".format(best_acc))
            print("Closing TF Session")

    def test(self, x_te, y_te):
        """Test routine"""

        with tf.Session() as sess:
            # Load the best model
            latest_checkpoint = tf.train.latest_checkpoint(
                self.save_dir)

            if tf.train.latest_checkpoint(self.save_dir) is not None:
                print("Restoring from {}...".format(
                    self.save_dir))
                self.saver_best.restore(
                    sess,
                    latest_checkpoint
                )

            # Test on the test data
            self.training = False
            num_test_b = len(x_te) // self.test_batch_size
            acc_list = []
            for idx_b in range(num_test_b):
                res = sess.run(
                    fetches={
                        "acc": self.acc,
                    },
                    feed_dict={
                        self.x_in: x_va[idx_b * self.test_batch_size: (idx_b + 1) * self.test_batch_size],
                        self.y_in: y_va[idx_b * self.test_batch_size: (idx_b + 1) * self.test_batch_size],
                    }
                )
                acc_list += [res['acc']]

            avg_acc = np.mean(acc_list)

            # Report (print) test result
            print("Test accuracy with the best model is {}".format(res_acc))

def main():

    print("Reading training data...")
    x_trva, y_trva = load_data("./data.h5", 12)

    # Randomly shuffle data and labels. IMPORANT: make sure the data and label
    # is shuffled with the same random indices so that they don't get mixed up!
    idx_shuffle = np.random.permutation(len(x_trva))
    x_trva = x_trva[idx_shuffle]
    y_trva = y_trva[idx_shuffle]

    # Change type to float32 and int64 since we are going to use that for
    # TensorFlow.
    x_trva = x_trva.astype("float32")
    y_trva = y_trva.astype("int64")

    # ----------------------------------------
    # Simply select the last 20% of the training data as validation dataset.
    num_tr = int(len(x_trva) * 0.8)
    x_tr = x_trva[:num_tr]
    x_va = x_trva[num_tr:]
    y_tr = y_trva[:num_tr]
    y_va = y_trva[num_tr:]
    frrn_a = FRRN((10, 720, 960, 3), 'A')

    # ----------------------------------------
    # Train
    # Run training
    frrn_a.train(x_tr, y_tr, x_va, y_va)

if __name__ == "__main__":

    main()

max_stars_repo_name: superpenshine/Full-Resolution-ResNet

max_issues_repo_name: superpenshine/Full-Resolution-ResNet

max_forks_repo_name: superpenshine/Full-Resolution-ResNet


text:
# Copyright 2019, FBPIC contributors
# Authors: Remi Lehe, Manuel Kirchen
# License: 3-Clause-BSD-LBNL
"""
This file is part of the Fourier-Bessel Particle-In-Cell code (FB-PIC)
It defines the charge and mass of standard particle types
"""
# Import constants from scipy
from scipy.constants import e, m_e, m_p, m_n, physical_constants

# Create dictionaries with species_type defined in openPMD 2
particle_charge = {
    'electron': -e,
    'positron': e,
    'proton': e,
    'anti-proton': -e,
    'neutron': 0,
    'anti-neutron': 0,
}
particle_mass = {
    'electron': m_e,
    'positron': m_e,
    'proton': m_p,
    'anti-proton': m_p,
    'neutron': m_n,
    'anti-neutron': m_n,
}

# Get mass of each element from periodictable
import periodictable
m_u = physical_constants['atomic mass constant'][0]
for el in periodictable.elements:
    particle_mass[ el.symbol ] = el.mass * m_u
    particle_charge[ el.symbol ] = 0

max_issues_repo_name: RemiLehe/fbpic

max_forks_repo_name: RemiLehe/fbpic

reason: from scipy


text:
import matplotlib
matplotlib.use('Agg')
import os
import numpy as np
import json
import caffe
from caffe import layers as L
from caffe import params as P
from vqa_data_provider_layer import VQADataProvider
from visualize_tools import drawgraph, exec_validation
import config

def qlstm(mode, batchsize, max_words_in_question, question_vocab_size):
    n = caffe.NetSpec()
    mode_str = json.dumps({'mode':mode, 'batchsize':batchsize})
    n.data, n.cont, n.img_feature, n.label = L.Python( \
        module='vqa_data_provider_layer', layer='VQADataProviderLayer', \
        param_str=mode_str, ntop=4 )

    n.embed_ba = L.Embed(n.data, input_dim=question_vocab_size, num_output=300, \
        weight_filler=dict(type='uniform',min=-0.08,max=0.08))
    n.embed = L.TanH(n.embed_ba)

    # LSTM1
    n.lstm1 = L.LSTM(\
                   n.embed, n.cont,\
                   recurrent_param=dict(\
                       num_output=1024,\
                       weight_filler=dict(type='uniform',min=-0.08,max=0.08),\
                       bias_filler=dict(type='constant',value=0)))
    tops1 = L.Slice(n.lstm1, ntop=max_words_in_question, slice_param={'axis':0})
    for i in xrange(max_words_in_question-1):
        n.__setattr__('slice_first'+str(i), tops1[int(i)])
        n.__setattr__('silence_data_first'+str(i), L.Silence(tops1[int(i)],ntop=0))
    n.lstm1_out = tops1[max_words_in_question-1]
    n.lstm1_reshaped = L.Reshape(n.lstm1_out,\
                          reshape_param=dict(\
                              shape=dict(dim=[-1,1024])))
    n.lstm1_reshaped_droped = L.Dropout(n.lstm1_reshaped,dropout_param={'dropout_ratio':0.3})
    n.lstm1_droped = L.Dropout(n.lstm1,dropout_param={'dropout_ratio':0.3})

    # LSTM2
    n.lstm2 = L.LSTM(\
                   n.lstm1_droped, n.cont,\
                   recurrent_param=dict(\
                       num_output=1024,\
                       weight_filler=dict(type='uniform',min=-0.08,max=0.08),\
                       bias_filler=dict(type='constant',value=0)))
    tops2 = L.Slice(n.lstm2, ntop=max_words_in_question, slice_param={'axis':0})
    for i in xrange(max_words_in_question-1):
        n.__setattr__('slice_second'+str(i), tops2[int(i)])
        n.__setattr__('silence_data_second'+str(i), L.Silence(tops2[int(i)],ntop=0))
    n.lstm2_out = tops2[max_words_in_question-1]
    n.lstm2_reshaped = L.Reshape(n.lstm2_out,\
                          reshape_param=dict(\
                              shape=dict(dim=[-1,1024])))
    n.lstm2_reshaped_droped = L.Dropout(n.lstm2_reshaped,dropout_param={'dropout_ratio':0.3})
    concat_botom = [n.lstm1_reshaped_droped, n.lstm2_reshaped_droped]
    n.lstm_12 = L.Concat(*concat_botom)

    n.q_emb_tanh_droped_resh = L.Reshape(n.lstm_12, \
        reshape_param=dict(shape=dict(dim=[-1,2048,1,1])))
    n.i_emb_tanh_droped_resh = L.Reshape(n.img_feature, \
        reshape_param=dict(shape=dict(dim=[-1,2048,1,1])))
    n.blcf = L.CompactBilinear(n.q_emb_tanh_droped_resh, n.i_emb_tanh_droped_resh, \
        compact_bilinear_param=dict(num_output=16000,sum_pool=False))
    n.blcf_sign_sqrt = L.SignedSqrt(n.blcf)
    n.blcf_sign_sqrt_l2 = L.L2Normalize(n.blcf_sign_sqrt)

    n.blcf_droped = L.Dropout(n.blcf_sign_sqrt_l2,dropout_param={'dropout_ratio':0.1})
    n.blcf_droped_resh = L.Reshape(n.blcf_droped,reshape_param=dict(shape=dict(dim=[-1,16000])))

    n.prediction = L.InnerProduct(n.blcf_droped_resh, num_output=config.NUM_OUTPUT_UNITS, \
        weight_filler=dict(type='xavier'))
    n.loss = L.SoftmaxWithLoss(n.prediction, n.label)
    return n.to_proto()

def make_answer_vocab(adic, vocab_size):
    """
    Returns a dictionary that maps words to indices.
    """
    adict = {'':0}
    nadict = {'':1000000}
    vid = 1
    for qid in adic.keys():
        answer_obj = adic[qid]
        answer_list = [ans['answer'] for ans in answer_obj]

        for q_ans in answer_list:
            # create dict
            if adict.has_key(q_ans):
                nadict[q_ans] += 1
            else:
                nadict[q_ans] = 1
                adict[q_ans] = vid
                vid +=1

    # debug
    nalist = []
    for k,v in sorted(nadict.items(), key=lambda x:x[1]):
        nalist.append((k,v))

    # remove words that appear less than once
    n_del_ans = 0
    n_valid_ans = 0
    adict_nid = {}
    for i, w in enumerate(nalist[:-vocab_size]):
        del adict[w[0]]
        n_del_ans += w[1]
    for i, w in enumerate(nalist[-vocab_size:]):
        n_valid_ans += w[1]
        adict_nid[w[0]] = i

    return adict_nid

def make_question_vocab(qdic):
    """
    Returns a dictionary that maps words to indices.
    """
    vdict = {'':0}
    vid = 1
    for qid in qdic.keys():
        # sequence to list
        q_str = qdic[qid]['qstr']
        q_list = VQADataProvider.seq_to_list(q_str)

        # create dict
        for w in q_list:
            if not vdict.has_key(w):
                vdict[w] = vid
                vid +=1

    return vdict

def make_vocab_files():
    """
    Produce the question and answer vocabulary files.
    """
    print 'making question vocab...', config.QUESTION_VOCAB_SPACE
    qdic, _ = VQADataProvider.load_data(config.QUESTION_VOCAB_SPACE)
    question_vocab = make_question_vocab(qdic)
    print 'making answer vocab...', config.ANSWER_VOCAB_SPACE
    _, adic = VQADataProvider.load_data(config.ANSWER_VOCAB_SPACE)
    answer_vocab = make_answer_vocab(adic, config.NUM_OUTPUT_UNITS)
    return question_vocab, answer_vocab

def main():
    if not os.path.exists('./result'):
        os.makedirs('./result')

    question_vocab, answer_vocab = {}, {}
    if os.path.exists('./result/vdict.json') and os.path.exists('./result/adict.json'):
        print 'restoring vocab'
        with open('./result/vdict.json','r') as f:
            question_vocab = json.load(f)
        with open('./result/adict.json','r') as f:
            answer_vocab = json.load(f)
    else:
        question_vocab, answer_vocab = make_vocab_files()
        with open('./result/vdict.json','w') as f:
            json.dump(question_vocab, f)
        with open('./result/adict.json','w') as f:
            json.dump(answer_vocab, f)

    print 'question vocab size:', len(question_vocab)
    print 'answer vocab size:', len(answer_vocab)

    with open('./result/proto_train.prototxt', 'w') as f:
        f.write(str(qlstm(config.TRAIN_DATA_SPLITS, config.BATCH_SIZE, \
            config.MAX_WORDS_IN_QUESTION, len(question_vocab))))

    with open('./result/proto_test.prototxt', 'w') as f:
        f.write(str(qlstm('val', config.BATCH_SIZE, \
            config.MAX_WORDS_IN_QUESTION, len(question_vocab))))

    caffe.set_device(config.GPU_ID)
    caffe.set_mode_gpu()
    solver = caffe.get_solver('./qlstm_solver.prototxt')

    train_loss = np.zeros(config.MAX_ITERATIONS)
    results = []

    for it in range(config.MAX_ITERATIONS):
        solver.step(1)

        # store the train loss
        train_loss[it] = solver.net.blobs['loss'].data

        if it % config.PRINT_INTERVAL == 0:
            print 'Iteration:', it
            c_mean_loss = train_loss[it-config.PRINT_INTERVAL:it].mean()
            print 'Train loss:', c_mean_loss
            if it % config.VALIDATE_INTERVAL == 0:
                print 'Validating...'
                solver.test_nets[0].save('./result/tmp.caffemodel')
                test_loss, acc_overall, acc_per_ques, acc_per_ans = exec_validation(config.GPU_ID, 'val', it=it)
                print 'Test loss:', test_loss
                print 'Accuracy:', acc_overall
                results.append([it, c_mean_loss, test_loss, acc_overall, acc_per_ques, acc_per_ans])
                best_result_idx = np.array([x[3] for x in results]).argmax()
                print 'Best accuracy of', results[best_result_idx][3], 'was at iteration', results[best_result_idx][0]
            drawgraph(results)

if __name__ == '__main__':
    main()

max_issues_repo_name: afcarl/vqa-mcb

max_forks_repo_name: afcarl/vqa-mcb


text:
import numpy as np

from ..item import Item
from .proto.step_pb2 import DType, ItemProto, Shape

def convert_dtype_to_proto(dtype: np.dtype) -> DType:
    if dtype == np.uint8:
        return DType.UINT8
    if dtype == np.int32:
        return DType.INT32
    if dtype == np.float32:
        return DType.FLOAT32
    raise ValueError(f"invalid dtype: {dtype}")

def convert_proto_to_dtype(dtype: DType) -> np.dtype:
    if dtype == DType.UINT8:
        return np.uint8
    if dtype == DType.INT32:
        return np.int32
    if dtype == DType.FLOAT32:
        return np.float32
    raise ValueError(f"invalid dtype: {dtype}")

def convert_item_to_proto(item: Item) -> ItemProto:
    if isinstance(item, (np.ndarray)):
        item = [item]
    elif isinstance(item, (int, float)):
        item = [np.array([item], dtype=np.float32)]
    assert isinstance(item, (list, tuple))

    shapes = [Shape(dim=el.shape) for el in item]
    data = [el.tobytes() for el in item]
    dtypes = [convert_dtype_to_proto(el.dtype) for el in item]

    return ItemProto(length=len(item), shape=shapes, data=data, dtype=dtypes)

def convert_proto_to_item(proto: ItemProto) -> Item:
    item = []
    for i in range(proto.length):
        dtype = convert_proto_to_dtype(proto.dtype[i])
        array = np.frombuffer(proto.data[i], dtype=dtype)
        item.append(np.reshape(array, proto.shape[i].dim))
    if len(item) == 1:
        return item[0]
    return item

max_stars_repo_name: takuseno/kiox

max_issues_repo_name: takuseno/kiox

max_forks_repo_name: takuseno/kiox


text:
# Copyright 2022 Dan Lim
#  Apache 2.0  (http://www.apache.org/licenses/LICENSE-2.0)

import numpy as np
import torch
import torch.nn as nn
import torch.nn.functional as F

from numba import jit

class AlignmentModule(nn.Module):
    """Alignment Learning Framework proposed for parallel TTS models in:

    https://arxiv.org/abs/2108.10447

    """

    def __init__(self, adim, odim):
        super().__init__()
        self.t_conv1 = nn.Conv1d(adim, adim, kernel_size=3, padding=1)
        self.t_conv2 = nn.Conv1d(adim, adim, kernel_size=1, padding=0)

        self.f_conv1 = nn.Conv1d(odim, adim, kernel_size=3, padding=1)
        self.f_conv2 = nn.Conv1d(adim, adim, kernel_size=3, padding=1)
        self.f_conv3 = nn.Conv1d(adim, adim, kernel_size=1, padding=0)

    def forward(self, text, feats, x_masks=None):
        """Calculate alignment loss.

        Args:
            text (Tensor): Batched text embedding (B, T_text, adim).
            feats (Tensor): Batched acoustic feature (B, T_feats, odim).
            x_masks (Tensor): Mask tensor (B, T_text).

        Returns:
            Tensor: Log probability of attention matrix (B, T_feats, T_text).

        """
        text = text.transpose(1, 2)
        text = F.relu(self.t_conv1(text))
        text = self.t_conv2(text)
        text = text.transpose(1, 2)

        feats = feats.transpose(1, 2)
        feats = F.relu(self.f_conv1(feats))
        feats = F.relu(self.f_conv2(feats))
        feats = self.f_conv3(feats)
        feats = feats.transpose(1, 2)

        dist = feats.unsqueeze(2) - text.unsqueeze(1)
        dist = torch.norm(dist, p=2, dim=3)
        score = -dist

        if x_masks is not None:
            x_masks = x_masks.unsqueeze(-2)
            score = score.masked_fill(x_masks, -np.inf)

        log_p_attn = F.log_softmax(score, dim=-1)
        return log_p_attn

@jit(nopython=True)
def _monotonic_alignment_search(log_p_attn):
    # https://arxiv.org/abs/2005.11129
    T_mel = log_p_attn.shape[0]
    T_inp = log_p_attn.shape[1]
    Q = np.full((T_inp, T_mel), fill_value=-np.inf)

    log_prob = log_p_attn.transpose(1, 0)  # -> (T_inp,T_mel)
    # 1.  Q <- init first row for all j
    for j in range(T_mel):
        Q[0, j] = log_prob[0, : j + 1].sum()

    # 2.
    for j in range(1, T_mel):
        for i in range(1, min(j + 1, T_inp)):
            Q[i, j] = max(Q[i - 1, j - 1], Q[i, j - 1]) + log_prob[i, j]

    # 3.
    A = np.full((T_mel,), fill_value=T_inp - 1)
    for j in range(T_mel - 2, -1, -1):  # T_mel-2, ..., 0
        # 'i' in {A[j+1]-1, A[j+1]}
        i_a = A[j + 1] - 1
        i_b = A[j + 1]
        if i_b == 0:
            argmax_i = 0
        elif Q[i_a, j] >= Q[i_b, j]:
            argmax_i = i_a
        else:
            argmax_i = i_b
        A[j] = argmax_i
    return A

def viterbi_decode(log_p_attn, text_lengths, feats_lengths):
    """Extract duration from an attention probability matrix

    Args:
        log_p_attn (Tensor): Batched log probability of attention
            matrix (B, T_feats, T_text).
        text_lengths (Tensor): Text length tensor (B,).
        feats_legnths (Tensor): Feature length tensor (B,).

    Returns:
        Tensor: Batched token duration extracted from `log_p_attn` (B, T_text).
        Tensor: Binarization loss tensor ().

    """
    B = log_p_attn.size(0)
    T_text = log_p_attn.size(2)
    device = log_p_attn.device

    bin_loss = 0
    ds = torch.zeros((B, T_text), device=device)
    for b in range(B):
        cur_log_p_attn = log_p_attn[b, : feats_lengths[b], : text_lengths[b]]
        viterbi = _monotonic_alignment_search(cur_log_p_attn.detach().cpu().numpy())
        _ds = np.bincount(viterbi)
        ds[b, : len(_ds)] = torch.from_numpy(_ds).to(device)

        t_idx = torch.arange(feats_lengths[b])
        bin_loss = bin_loss - cur_log_p_attn[t_idx, viterbi].mean()
    bin_loss = bin_loss / B
    return ds, bin_loss

@jit(nopython=True)
def _average_by_duration(ds, xs, text_lengths, feats_lengths):
    B = ds.shape[0]
    xs_avg = np.zeros_like(ds)
    ds = ds.astype(np.int32)
    for b in range(B):
        t_text = text_lengths[b]
        t_feats = feats_lengths[b]
        d = ds[b, :t_text]
        d_cumsum = d.cumsum()
        d_cumsum = [0] + list(d_cumsum)
        x = xs[b, :t_feats]
        for n, (start, end) in enumerate(zip(d_cumsum[:-1], d_cumsum[1:])):
            if len(x[start:end]) != 0:
                xs_avg[b, n] = x[start:end].mean()
            else:
                xs_avg[b, n] = 0
    return xs_avg

def average_by_duration(ds, xs, text_lengths, feats_lengths):
    """Average frame-level features into token-level according to durations

    Args:
        ds (Tensor): Batched token duration (B, T_text).
        xs (Tensor): Batched feature sequences to be averaged (B, T_feats).
        text_lengths (Tensor): Text length tensor (B,).
        feats_lengths (Tensor): Feature length tensor (B,).

    Returns:
        Tensor: Batched feature averaged according to the token duration (B, T_text).

    """
    device = ds.device
    args = [ds, xs, text_lengths, feats_lengths]
    args = [arg.detach().cpu().numpy() for arg in args]
    xs_avg = _average_by_duration(*args)
    xs_avg = torch.from_numpy(xs_avg).to(device)
    return xs_avg

max_stars_repo_name: imdanboy/espnet

max_issues_repo_name: imdanboy/espnet

max_forks_repo_name: imdanboy/espnet

reason: import numpy,from numba