text:
#   Copyright (c) 2019 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 copy
import numpy as np

from .framework import Variable, default_main_program, default_startup_program, _in_dygraph_mode, _current_expected_place
from . import unique_name
from .param_attr import ParamAttr, WeightNormParamAttr
from . import core

class LayerHelperBase(object):
    def __init__(self, name, layer_type):
        self._layer_type = layer_type
        self._name = name

    @property
    def name(self):
        return self._name

    @property
    def layer_type(self):
        return self._layer_type

    @property
    def main_program(self):
        return default_main_program()

    @property
    def startup_program(self):
        return default_startup_program()

    def to_variable(self, value, block=None):
        """convert value to variable

            Args:
                value: value to be convert
                block: the block of the variable

        Return Variable construct from value
        """
        if isinstance(value, np.ndarray):
            assert _in_dygraph_mode(
            ), "to_variable could only be called in dygraph mode"

            if not block:
                block = default_main_program().current_block()
            py_var = Variable(
                block,
                type=core.VarDesc.VarType.LOD_TENSOR,
                name=None,
                shape=value.shape,
                dtype=value.dtype)
            var = py_var._ivar.value()
            tensor = var.get_tensor()
            tensor.set(value, _current_expected_place())
            return py_var
        elif isinstance(value, Variable):
            return value

    def _create_weight_normalize(self, attr, shape, dtype):
        from .layers import elementwise_mul, elementwise_div, reshape

        # Remove these ops when LayerHelper and layers support indicating
        # program and block.
        def __norm_op(x,
                      out=None,
                      p=2,
                      dim=None,
                      keep_dim=False,
                      block=self.startup_program.global_block()):
            if out is None:
                out = block.create_var(
                    name=unique_name.generate(".".join(
                        [self.name, 'weight_norm_norm'])),
                    dtype=dtype,
                    persistable=False)
            abs_out = block.create_var(
                name=unique_name.generate(".".join(
                    [self.name, 'weight_norm_abs'])),
                dtype=dtype,
                persistable=False)
            block.append_op(
                type='abs', inputs={'X': x}, outputs={'Out': abs_out})
            pow_out = block.create_var(
                name=unique_name.generate(".".join(
                    [self.name, 'weight_norm_pow'])),
                dtype=dtype,
                persistable=False)
            block.append_op(
                type='pow',
                inputs={'X': abs_out},
                outputs={'Out': pow_out},
                attrs={'factor': float(p)})
            sum_out = block.create_var(
                name=unique_name.generate(".".join(
                    [self.name, 'weight_norm_sum'])),
                dtype=dtype,
                persistable=False)
            block.append_op(
                type='reduce_sum',
                inputs={'X': pow_out},
                outputs={'Out': sum_out},
                attrs={
                    'dim': dim,
                    'keep_dim': keep_dim,
                    'reduce_all': True if dim is None else False
                })
            block.append_op(
                type='pow',
                inputs={'X': sum_out},
                outputs={'Out': out},
                attrs={'factor': 1. / p})
            return out

        def __reshape_op(x,
                         shape,
                         out=None,
                         block=self.startup_program.global_block()):
            if out is None:
                out = block.create_var(
                    name=unique_name.generate(".".join(
                        [self.name, 'weight_norm_reshape'])),
                    dtype=dtype,
                    persistable=False)
            block.append_op(
                type='reshape',
                inputs={'X': x},
                outputs={'Out': out},
                attrs={'shape': shape})
            return out

        def __transpose_op(x,
                           axis,
                           out=None,
                           block=self.startup_program.global_block()):
            if out is None:
                out = block.create_var(
                    name=unique_name.generate(".".join(
                        [self.name, 'weight_norm_transpose'])),
                    dtype=dtype,
                    persistable=False)
            block.append_op(
                type='transpose',
                inputs={'X': x},
                outputs={'Out': out},
                attrs={'axis': axis})
            return out

        def __norm_except_dim(x,
                              out=None,
                              dim=None,
                              block=self.startup_program.global_block()):
            """Computes the norm over all dimensions except dim"""
            if out is None:
                out = block.create_var(
                    name=unique_name.generate(".".join(
                        [self.name, 'weight_norm_norm'])),
                    dtype=dtype,
                    persistable=False)
            if dim is None:
                __norm_op(x, out, dim=dim, block=block)
            elif dim == 0:
                out_shape = [x.shape[0]] + [1] * (len(x.shape) - 1)
                reshape = __reshape_op(x, shape=[x.shape[0], -1], block=block)
                norm = __norm_op(reshape, dim=1, block=block)
                __reshape_op(norm, out=out, shape=out_shape, block=block)
            elif dim == len(x.shape) - 1:
                out_shape = [1] * (len(x.shape) - 1) + [x.shape[-1]]
                reshape = __reshape_op(x, shape=[-1, x.shape[-1]], block=block)
                norm = __norm_op(reshape, dim=0, block=block)
                __reshape_op(norm, out=out, shape=out_shape, block=block)
            else:
                perm = list(range(len(x.shape)))
                perm[0], perm[dim] = dim, 0
                transpose = __transpose_op(x, perm, block=block)
                norm = __norm_op(transpose, dim=0, block=block)
                __transpose_op(norm, perm, out=out, block=block)
            return out

        def __weight_normalize(g, v, dim):
            """Calculations for weight normalization"""
            norm = __norm_except_dim(
                v, dim=dim, block=self.main_program.current_block())
            scale = elementwise_div(
                x=g, y=norm)  # The shapes of g and norm are the same.
            # Currently, elementwise_mul only support broadcast when the shape
            # of y is a subset of the shape of x. Thus, we reshape y to squeeze
            # to achive the subset.
            w = elementwise_mul(
                x=v,
                y=scale if dim is None else reshape(
                    x=scale, shape=[v.shape[dim]]),
                axis=-1 if dim is None else dim)
            # To serialize the original parameter for inference, maybe a
            # parameter rather than a variable should be returned.
            return w

        g_param_attr = copy.deepcopy(attr)
        g_param_attr.name = attr.name + '_g'
        g_param_shape = [1] * len(shape)
        if attr.dim is not None:
            g_param_shape[attr.dim] = shape[attr.dim]
        v_param_attr = copy.deepcopy(attr)
        v_param_attr.name = attr.name + '_v'
        v_param_shape = shape

        # Add to startup_program to initialize g and v.
        # Try to reconstruct the initializer of w by initializing g and v.
        # Set the initializers of g and v as below, then the distribution
        # of w is the same as initializing w with the given initializer.
        # For Data-Dependent Initialization, please compute the init-values
        # of g and v in external and then feed the values to g and v by
        # executing an extra program.
        g_param = self.startup_program.global_block().create_parameter(
            dtype=dtype,
            shape=g_param_shape,
            **g_param_attr._to_kwargs(with_initializer=False))
        v_param = self.startup_program.global_block().create_parameter(
            dtype=dtype,
            shape=v_param_shape,
            **v_param_attr._to_kwargs(with_initializer=True))
        __norm_except_dim(
            x=v_param,
            out=g_param,
            dim=attr.dim,
            block=self.startup_program.global_block())

        # Add weight normalization to main_program
        g_param = self.main_program.global_block().create_parameter(
            dtype=dtype, shape=g_param_shape, **g_param_attr._to_kwargs())
        v_param = self.main_program.global_block().create_parameter(
            dtype=dtype, shape=v_param_shape, **v_param_attr._to_kwargs())
        w_param = __weight_normalize(g_param, v_param, dim=attr.dim)
        return w_param

    # TODO: hide the func after we move the layers to Layers
    def create_parameter(self,
                         attr,
                         shape,
                         dtype,
                         is_bias=False,
                         default_initializer=None):
        """Create parameters for this layers.

           Args:
               attr: [ParamAttr] should be the parameter attribute for this parameter
               shape: shape of the paramter
               dtype: data type of this parameter
               is_bias: if this is a bias parameter
               default_initializer: set the default initializer for this parameter

        Returns created parameter Variable.
        """
        # Deepcopy the attr so that parameters can be shared in program
        attr = copy.deepcopy(attr)
        attr = ParamAttr._to_attr(attr)
        if not attr:
            return None
        assert isinstance(attr, ParamAttr)
        suffix = 'b' if is_bias else 'w'
        if attr.name is None:
            attr.name = unique_name.generate(".".join([self.name, suffix]))

        if default_initializer is None and attr.initializer is None:
            if isinstance(dtype, core.VarDesc.VarType):
                if dtype != core.VarDesc.VarType.FP32 and \
                        dtype != core.VarDesc.VarType.FP64 and \
                        dtype != core.VarDesc.VarType.FP16:
                    raise TypeError(
                        "Can not create parameter with default initializer when dtype is not float type. Set default_initializer to fit the parameter dtype!"
                    )
            else:
                if not (dtype.startswith("float") or dtype == "double"):
                    raise TypeError(
                        "Can not create parameter with default initializer when dtype is not float type. Set default_initializer to fit the parameter dtype!"
                    )
            if is_bias:
                attr._set_default_bias_initializer()
            else:
                attr._set_default_param_initializer()
        else:
            attr._set_default_initializer(default_initializer)

        # If weight normalization is set, insert extra parameters and ops.
        # Refer to https://arxiv.org/pdf/1602.07868.pdf
        if isinstance(attr, WeightNormParamAttr):
            param = self._create_weight_normalize(attr, shape, dtype)
            WeightNormParamAttr.params_with_weight_norm.append(param)
            return param
        if _in_dygraph_mode():
            # In dygraph mode, we want the returned parameter to be
            # initialized so that it can be used imperatively.
            return self.main_program.global_block().create_parameter(
                dtype=dtype,
                shape=shape,
                **attr._to_kwargs(with_initializer=True))
        else:
            self.startup_program.global_block().create_parameter(
                dtype=dtype,
                shape=shape,
                **attr._to_kwargs(with_initializer=True))
            return self.main_program.global_block().create_parameter(
                dtype=dtype, shape=shape, **attr._to_kwargs())

    def create_variable_for_type_inference(self, dtype, stop_gradient=False):
        """Create a temporary variable that should be type inferred layer.

        Note:
            The default type will be set to LOD_TENSOR. However, when
            the var is used as operator output, its type will be updated
            based on operator's `VarTypeInference` implementation in
            infer_var_type.
        """
        return self.main_program.current_block().create_var(
            name=unique_name.generate(".".join([self.name, 'tmp'])),
            dtype=dtype,
            type=core.VarDesc.VarType.LOD_TENSOR,
            persistable=False,
            stop_gradient=stop_gradient)

    def create_variable(self, *args, **kwargs):
        """Create Variable for this layers.
        Returns created Variable.
        """
        return self.main_program.current_block().create_var(*args, **kwargs)

    def create_global_variable(self, persistable=False, *args, **kwargs):
        """
        create global variable, note that there is no initializer for this global variable.
        Args:
            persistable(bool): True if it is a checkpoint value.
            *args: See create_var's documentation
            **kwargs: See create_var's documentation

        Returns(Variable): the created variable.
        """
        return self.main_program.global_block().create_var(
            *args, persistable=persistable, **kwargs)

    def create_or_get_global_variable(self, name, *args, **kwargs):
        """
        Creates a global variable if not exists and returns the variable and
        a boolean flag which is true when it is a new variable.
        """
        if self.main_program.global_block().has_var(name):
            return self.main_program.global_block().var(name), False
        else:
            return self.create_global_variable(name=name, *args, **kwargs), True

    def set_variable_initializer(self, var, initializer):
        """Set target Variable's initializer

           Args:
               var: target Variable
               initializer: initializer to use
        """
        assert isinstance(var, Variable)
        if _in_dygraph_mode():
            initializer(var, var.block)
        else:
            self.startup_program.global_block().create_var(
                name=var.name,
                type=var.type,
                dtype=var.dtype,
                shape=var.shape,
                persistable=True,
                initializer=initializer)


text:
import numpy as np
import time
import torch
import torch.nn as nn
import torch.autograd
import h5py
import torch.optim as optim
import scipy.io
from torch.autograd import Variable
import torch.optim as optim
from enum import Enum

from HeatEquation.Dataset.Baseline import HeatEquationDataset
from Schrodinger.Schrodinger2D_DeepHPM import SchrodingerNet

import matplotlib.pyplot as plt
import torch.utils.data.distributed
import horovod.torch as hvd
from tensorboardX import SummaryWriter
from argparse import ArgumentParser
import os
import sys
import pathlib
import torch.nn.functional as F

class HeatEquationNet(SchrodingerNet):
    def __init__(self, numLayers, numFeatures, numLayers_hpm, numFeatures_hpm, lb, ub, samplingX, samplingY, activation=torch.tanh, activation_hpm=F.relu):
        """
        This function creates the components of the Neural Network and saves the datasets
        :param x0: Position x at time zero
        :param u0: Real Part of the solution at time 0 at position x
        :param v0: Imaginary Part of the solution at time 0 at position x
        :param tb: Time Boundary
        :param X_f: Training Data for partial differential equation
        :param layers: Describes the structure of Neural Network
        :param lb: Value of the lower bound in space
        :param ub: Value of the upper bound in space
        """
        torch.manual_seed(1234)
        super(HeatEquationNet, self).__init__(numLayers, numFeatures, numLayers_hpm, numFeatures_hpm, lb, ub, samplingX, samplingY, activation, activation_hpm)

    def init_layers(self):
        """
        This function creates the torch layers and initialize them with xavier
        :param self:
        :return:
        """

        self.lin_layers.append(nn.Linear(3, self.numFeatures))
        for _ in range(self.numLayers):
            inFeatures = self.numFeatures
            self.lin_layers.append(nn.Linear(inFeatures, self.numFeatures))
        self.lin_layers.append(nn.Linear(inFeatures, 1))

        for m in self.lin_layers:
            if isinstance(m, nn.Linear):
                nn.init.xavier_normal_(m.weight)
                nn.init.constant_(m.bias, 0)

        self.lin_layers_hpm.append(nn.Linear(6, self.numFeatures))
        for _ in range(self.numLayers_hpm):
            inFeatures = self.numFeatures_hpm
            self.lin_layers_hpm.append(nn.Linear(inFeatures, self.numFeatures_hpm))
        self.lin_layers_hpm.append(nn.Linear(inFeatures, 1))

        for m in self.lin_layers_hpm:
            if isinstance(m, nn.Linear):
                nn.init.xavier_normal_(m.weight)
                nn.init.constant_(m.bias, 0)

    def net_uv(self, x, y, t):
        """
        Function that calculates the nn output at postion x,y at time t
        :param x: position
        :param y: position
        :param t: time
        :return: Solutions and their gradients
        """

        dim = x.shape[0]
        x = Variable(x, requires_grad=True).cuda()
        y = Variable(y, requires_grad=True).cuda()
        t = Variable(t, requires_grad=True).cuda()

        X = torch.stack([x, y, t], 1)
        u = self.forward(X)

        grads = torch.ones([dim]).cuda()

        # huge change to the tensorflow implementation this function returns all neccessary gradients
        u_x = torch.autograd.grad(u, x, create_graph=True, grad_outputs=grads)[0]
        u_y = torch.autograd.grad(u, y, create_graph=True, grad_outputs=grads)[0]
        u_x = u_x.reshape([dim])
        u_y = u_x.reshape([dim])
        u_t = torch.autograd.grad(u, t, create_graph=True, grad_outputs=grads)[0]
        u_xx = torch.autograd.grad(u_x, x, create_graph=True, grad_outputs=grads)[0]
        u_yy = torch.autograd.grad(u_y, y, create_graph=True, grad_outputs=grads)[0]

        u_t = u_t.reshape([dim])
        u_xx = u_xx.reshape([dim])
        u_yy = u_yy.reshape([dim])

        return u, u_yy, u_xx, u_t

    def net_pde(self, x, y, t, gamma=1.):
        """
        Calculates the quality of the pde estimation
        :param x postion x
        :param y postion y
        :param t time t
        """
        u, u_yy, u_xx, u_t = self.net_uv(x, y, t)
        x = x.view(-1)
        y = y.view(-1)

        X = torch.stack([x, y, t, u, u_yy, u_xx], 1)

        f = -1 * u_t - self.forward_hpm(X)

        return u, f

    def pde_loss(self, x0, y0, t0, u0, v0, xf, yf, tf):
        """
        Returns the quality of the net
        """

        x0 = x0.view(-1)
        y0 = y0.view(-1)
        t0 = t0.view(-1)
        xf = xf.view(-1)
        yf = yf.view(-1)
        tf = tf.view(-1)

        n0 = x0.shape[0]

        inputX = torch.cat([x0, xf])
        inputY = torch.cat([y0, yf])
        inputT = torch.cat([t0, tf])

        u, f = self.net_pde(inputX, inputY, inputT)

        solU = u[:n0]

        u0 = u0.view(-1)

        pdeLoss = torch.mean((solU - u0) ** 2) + torch.mean(f ** 2)
        return pdeLoss

    def hpm_loss(self, x, y, t, Ex_u):
        """
        Returns the quality HPM net
        """

        x = x.view(-1)
        y = y.view(-1)
        t = t.view(-1)

        u, f = self.net_pde(x, y, t)

        Ex_u = Ex_u.view(-1)

        return torch.mean((u - Ex_u) ** 2) + torch.mean(f ** 2)

def writeIntermediateState(timeStep, model, epoch, nx, ny, fileWriter,csystem):
    if fileWriter:
        x, y, t = HeatEquationDataset.getInput(timeStep,csystem)
        x = torch.Tensor(x).float().cuda()
        y = torch.Tensor(y).float().cuda()
        t = torch.Tensor(t).float().cuda()

        inputX = torch.stack([x, y, t], 1)
        u = model.forward(inputX).detach().cpu().numpy()

        fig = plt.figure()
        plt.imshow(u, cmap='jet')
        plt.colorbar()
        fileWriter.add_figure('u_' + str(timeStep), fig, epoch)
        plt.close(fig)

def valLoss(model, timeStep, csystem):
    x, y, t = HeatEquationDataset.getInput(timeStep,csystem)
    x = torch.Tensor(x).float().cuda()
    y = torch.Tensor(y).float().cuda()
    t = torch.Tensor(t).float().cuda()

    inputX = torch.stack([x, y, t], 1)
    uPred = model.forward(inputX).detach().cpu().numpy()

    # load label data
    uVal, vVal = HeatEquationDataset.loadFrame(pData, timeStep)
    uVal = np.array(uVal)

    return np.mean((uVal - uPred) ** 2)

def writeValidationLoss(model, writer, timeStep, epoch, csystem):
    if writer:
        loss = valLoss(model, timeStep,csystem)
        writer.add_scalar("ValidationLoss_" + str(timeStep), loss, epoch)

def save_checkpoint(model, optimizer, path, epoch):
    state = {
        'model': model.state_dict(),
        'optimizer': optimizer.state_dict(),
    }
    torch.save(state, path + 'model_' + str(epoch))

def load_checkpoint(model, optimizer, path):
    checkpoint = torch.load(path)
    model.load_state_dict(checkpoint['model'])
    optimizer.load_state_dict(checkpoint['optimizer'])

if __name__ == "__main__":

    # Initialize Horovod
    hvd.init()

    # Pin GPU to be used to process local rank (one GPU per process)
    torch.cuda.set_device(hvd.local_rank())

    # static parameter
    nx = 480
    ny = 640
    nt = 1000
    xmin = -2.5
    xmax = 2.5
    ymin = -2.5
    ymax = 2.5
    dt = 0.016 # roughly 60 Hz

    coordinateSystem = {"x_lb": xmin, "x_ub": xmax, "y_lb": ymin, "y_ub" : ymax, "nx": nx , "ny": ny, "nt": nt, "dt": dt}

    pData = 'data/UKD/2018_002/20180115_114053'

    parser = ArgumentParser()
    parser.add_argument("--identifier", dest="identifier", type=str)
    parser.add_argument("--batchsize", dest="batchsize", type=int)
    parser.add_argument("--numbatches", dest="numBatches", type=int)
    parser.add_argument("--initsize", dest="initSize", type=int)
    parser.add_argument("--numlayers", dest="numLayers", type=int)
    parser.add_argument("--numfeatures", dest="numFeatures", type=int)
    parser.add_argument("--numlayers_hpm", dest="numLayers_hpm", type=int)
    parser.add_argument("--numfeatures_hpm", dest="numFeatures_hpm", type=int)
    parser.add_argument("--epochsPDE", dest="epochsPDE", type=int)
    parser.add_argument("--alpha",dest="alpha",type=float)

    args = parser.parse_args()

    if hvd.rank() == 0:
        print("-" * 10 + "-" * len(args.identifier) + "-" * 10)
        print("-" * 10 +  args.identifier + "-" * 10)
        print("-" * 10 + "-" * len(args.identifier) + "-" * 10)

    print("Rank",hvd.rank(),"Local Rank", hvd.local_rank())

    #adapter of commandline parameters
    modelPath = 'runs/models/' + args.identifier + '/'
    logdir = 'runs/logs/' + args.identifier
    batchSizePDE = args.batchsize
    useGPU = True
    numBatches = args.numBatches
    initSize = args.initSize
    numLayers = args.numLayers
    numFeatures = args.numFeatures
    numLayers_hpm = args.numLayers_hpm
    numFeatures_hpm = args.numFeatures_hpm
    numEpochsPDE = args.epochsPDE
    activation = torch.tanh

    #create modelpath
    if hvd.rank() == 0:
        pathlib.Path(modelPath).mkdir(parents=True, exist_ok=True)
    # create logWriter
    log_writer = SummaryWriter(logdir) if hvd.rank() == 0 else None

    # create dataset
    ds = HeatEquationDataset(pData, coordinateSystem,  numBatches, batchSizePDE, useGPU=True)

    # Partition dataset among workers using DistributedSampler
    train_sampler = torch.utils.data.distributed.DistributedSampler(ds, num_replicas=hvd.size(), rank=hvd.rank())
    train_loader = torch.utils.data.DataLoader(ds, batch_size=1, sampler=train_sampler)

    model = HeatEquationNet(numLayers, numFeatures, numLayers_hpm, numFeatures_hpm, ds.lb, ds.ub, numOfEnergySamplingPointsX, numOfEnergySamplingPointsY, torch.tanh).cuda()

    optimizer = optim.Adam(model.parameters(), lr=3e-5)
    optimizer = hvd.DistributedOptimizer(optimizer,
                                         named_parameters=model.named_parameters(),
                                         backward_passes_per_step=1)

    for epoch in range(numEpochsPDE):

        for x, y, t, Ex_u in train_loader:
            optimizer.zero_grad()

            # calculate loss

            loss = model.hpm_loss(x,
                                  y,
                                  t,
                                  Ex_u)
            loss.backward()
            optimizer.step()

        if epoch % 30 == 0:
            writeIntermediateState(0, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(100, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(200, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(250, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(300, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(400, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeIntermediateState(500, model, epoch, nx, ny, log_writer,coordinateSystem)
            writeValidationLoss(model, log_writer, 250, epoch,coordinateSystem)
            writeValidationLoss(model, log_writer, 500, epoch,coordinateSystem)
            sys.stdout.flush()

            print("PDE Loss at Epoch: ", epoch + 1, loss.item())
            if log_writer:
                log_writer.add_histogram('First Layer Grads', model.lin_layers[0].weight.grad.view(-1, 1), epoch)
                save_checkpoint(model, optimizer, modelPath, epoch)

ext: py

max_stars_repo_name: maxxxzdn/PDElearning

max_issues_repo_name: maxxxzdn/PDElearning

max_forks_repo_name: maxxxzdn/PDElearning

reason: import numpy,import scipy


text:
import numpy as np

class LogisticRegressionModel(object):
    def __init__(self, weights, b, x_mean=None, x_scale=None, sta=None, phase=None):

        self.weights = weights
        self.b = b

        if x_mean is None:
            x_mean = np.zeros(weights.shape)
        self.x_mean = x_mean

        if x_scale is None:
            x_scale = np.ones(weights.shape)
        self.x_scale = x_scale

        self.sta = sta
        self.phase = phase

    def predict_prob(self, x):
        centered = (x-self.x_mean) / self.x_scale
        log_odds = np.dot(centered, self.weights) + self.b
        return 1.0 / (1.0 + np.exp(-log_odds))

max_stars_repo_name: davmre/sigvisa

max_issues_repo_name: davmre/sigvisa

max_forks_repo_name: davmre/sigvisa

reason: import numpy


text:
from typing import Optional, Tuple

import cv2
import numpy as np

from l5kit.geometry import transform_point, transform_points

PREDICTED_POINTS_COLOR = (0, 255, 255)
TARGET_POINTS_COLOR = (255, 0, 255)
REFERENCE_TRAJ_COLOR = (255, 255, 0)
#  Arrows represent position + orientation.
ARROW_LENGTH_IN_PIXELS = 2
ARROW_THICKNESS_IN_PIXELS = 1
ARROW_TIP_LENGTH_IN_PIXELS = 1.8

def draw_arrowed_line(on_image: np.ndarray, position: np.ndarray, yaw: float, rgb_color: Tuple[int, int, int]) -> None:
    """
    Draw a single arrowed line in an RGB image
    Args:
        on_image (np.ndarray): the RGB image to draw onto
        position (np.ndarray): the starting position of the arrow
        yaw (float): the arrow orientation
        rgb_color (Tuple[int, int, int]): the arrow color

    Returns: None

    """
    start_pixel = np.array(position[:2])

    rot = cv2.getRotationMatrix2D((0, 0), np.degrees(-yaw), 1.0)  # minus here because of cv2 rotations convention
    end_pixel = start_pixel + transform_point(np.asarray([ARROW_LENGTH_IN_PIXELS, 0]), rot)

    cv2.arrowedLine(
        on_image,
        tuple(start_pixel.astype(np.int32)),
        tuple(end_pixel.astype(np.int32)),
        rgb_color,
        thickness=ARROW_THICKNESS_IN_PIXELS,
        tipLength=ARROW_TIP_LENGTH_IN_PIXELS,
    )

def draw_trajectory(
    on_image: np.ndarray,
    positions: np.ndarray,
    rgb_color: Tuple[int, int, int],
    radius: int = 1,
    yaws: Optional[np.ndarray] = None,
) -> None:
    """
    Draw a trajectory on oriented arrow onto an RGB image
    Args:
        on_image (np.ndarray): the RGB image to draw onto
        positions (np.ndarray): pixel coordinates in the image space (not displacements) (Nx2)
        rgb_color (Tuple[int, int, int]): the trajectory RGB color
        radius (int): radius of the circle
        yaws (Optional[np.ndarray]): yaws in radians (N) or None to disable yaw visualisation

    Returns: None

    """
    if yaws is not None:
        assert len(yaws) == len(positions)
        for pos, yaw in zip(positions, yaws):
            pred_waypoint = pos[:2]
            pred_yaw = float(yaw[0])
            draw_arrowed_line(on_image, pred_waypoint, pred_yaw, rgb_color)
    else:
        for pos in positions:
            pred_waypoint = pos[:2]
            cv2.circle(on_image, tuple(pred_waypoint.astype(np.int)), radius, rgb_color, -1)

def draw_reference_trajectory(on_image: np.ndarray, world_to_pixel: np.ndarray, positions: np.ndarray) -> None:
    """
    Draw a trajectory (as points) onto the image
    Args:
        on_image (np.ndarray): the RGB image to draw onto
        world_to_pixel (np.ndarray): 3x3 matrix from meters to ego pixel space
        positions (np.ndarray): positions as 2D absolute meters coordinates

    Returns: None

    """
    positions_in_pixel_space = transform_points(np.array(positions), world_to_pixel)
    #  we clip the positions to be within the image
    mask = np.all(positions_in_pixel_space > (0.0, 0.0), 1) * np.all(positions_in_pixel_space < on_image.shape[:2], 1)
    positions_in_pixel_space = positions_in_pixel_space[mask]
    for pos in positions_in_pixel_space:
        cv2.circle(on_image, tuple(np.floor(pos).astype(np.int32)), 1, REFERENCE_TRAJ_COLOR, -1)


text:
"""
Unified methods of projection between sicd pixel coordinates,
some ortho-rectified pixel grid coordinates, and geophysical coordinates
"""

__classification__ = "UNCLASSIFIED"
__author__ = "Thomas McCullough"

import logging

import numpy
from sarpy.io.complex.sicd_elements.SICD import SICDType
from sarpy.geometry.geocoords import geodetic_to_ecf, ecf_to_geodetic, wgs_84_norm
from sarpy.processing.rational_polynomial import SarpyRatPolyError, \
    get_rational_poly_2d, get_rational_poly_3d, CombinedRationalPolynomial

logger = logging.getLogger(__name__)

_PIXEL_METHODOLOGY = ('MAX', 'MIN', 'MEAN', 'GEOM_MEAN')

class ProjectionHelper(object):
    """
    Abstract helper class which defines the projection interface for
    ortho-rectification usage for a sicd type object.
    """

    __slots__ = ('_sicd', '_row_spacing', '_col_spacing', '_default_pixel_method')

    def __init__(self, sicd, row_spacing=None, col_spacing=None, default_pixel_method='GEOM_MEAN'):
        r"""

        Parameters
        ----------
        sicd : SICDType
            The sicd object
        row_spacing : None|float
            The row pixel spacing. If not provided, this will default according
            to `default_pixel_method`.
        col_spacing : None|float
            The row pixel spacing. If not provided, this will default according
            to `default_pixel_method`.
        default_pixel_method : str
            Must be one of ('MAX', 'MIN', 'MEAN', 'GEOM_MEAN'). This determines
            the default behavior for row_spacing/col_spacing. The default value for
            row/column spacing will be the implied function applied to the range
            and azimuth ground resolution. Note that geometric mean is defined as
            :math:`\sqrt(x*x + y*y)`
        """

        self._row_spacing = None
        self._col_spacing = None
        default_pixel_method = default_pixel_method.upper()
        if default_pixel_method not in _PIXEL_METHODOLOGY:
            raise ValueError(
                'default_pixel_method got invalid value {}. Must be one '
                'of {}'.format(default_pixel_method, _PIXEL_METHODOLOGY))
        self._default_pixel_method = default_pixel_method

        if not isinstance(sicd, SICDType):
            raise TypeError('sicd must be a SICDType instance. Got type {}'.format(type(sicd)))
        if not sicd.can_project_coordinates():
            raise ValueError('Ortho-rectification requires the SICD ability to project coordinates.')
        sicd.define_coa_projection(overide=False)
        self._sicd = sicd
        self.row_spacing = row_spacing
        self.col_spacing = col_spacing

    @property
    def sicd(self):
        """
        SICDType: The sicd structure.
        """

        return self._sicd

    @property
    def row_spacing(self):
        """
        float: The row pixel spacing
        """

        return self._row_spacing

    @row_spacing.setter
    def row_spacing(self, value):
        """
        Set the row pixel spacing value. Setting to None will result in a
        default value derived from the SICD structure being used.

        Parameters
        ----------
        value : None|float

        Returns
        -------
        None
        """

        if value is None:
            if self.sicd.RadarCollection.Area is None:
                self._row_spacing = self._get_sicd_ground_pixel()
            else:
                self._row_spacing = self.sicd.RadarCollection.Area.Plane.XDir.LineSpacing
        else:
            value = float(value)
            if value <= 0:
                raise ValueError('row pixel spacing must be positive.')
            self._row_spacing = float(value)

    @property
    def col_spacing(self):
        """
        float: The column pixel spacing
        """

        return self._col_spacing

    @col_spacing.setter
    def col_spacing(self, value):
        """
        Set the col pixel spacing value. Setting to NOne will result in a
        default value derived from the SICD structure being used.

        Parameters
        ----------
        value : None|float

        Returns
        -------
        None
        """

        if value is None:
            if self.sicd.RadarCollection.Area is None:
                self._col_spacing = self._get_sicd_ground_pixel()
            else:
                self._col_spacing = self.sicd.RadarCollection.Area.Plane.YDir.SampleSpacing
        else:
            value = float(value)
            if value <= 0:
                raise ValueError('column pixel spacing must be positive.')
            self._col_spacing = float(value)

    def _get_sicd_ground_pixel(self):
        """
        Gets the SICD ground pixel size.

        Returns
        -------
        float
        """

        ground_row_ss, ground_col_ss = self.sicd.get_ground_resolution()
        if self._default_pixel_method == 'MIN':
            return min(ground_row_ss, ground_col_ss)
        elif self._default_pixel_method == 'MAX':
            return max(ground_row_ss, ground_col_ss)
        elif self._default_pixel_method == 'MEAN':
            return 0.5*(ground_row_ss + ground_col_ss)
        elif self._default_pixel_method == 'GEOM_MEAN':
            return float(numpy.sqrt(ground_row_ss*ground_row_ss + ground_col_ss*ground_col_ss))
        else:
            raise ValueError('Got unhandled default_pixel_method {}'.format(self._default_pixel_method))

    @staticmethod
    def _reshape(array, final_dimension):
        """
        Reshape the input so that the output is two-dimensional with final
        dimension given by `final_dimension`.

        Parameters
        ----------
        array : numpy.ndarray|list|tuple
        final_dimension : int

        Returns
        -------
        (numpy.ndarray, tuple)
            The reshaped data array and original shape.
        """

        if not isinstance(array, numpy.ndarray):
            array = numpy.array(array, dtype=numpy.float64)
        if array.ndim < 1 or array.shape[-1] != final_dimension:
            raise ValueError(
                'ortho_coords must be at least one dimensional with final dimension '
                'of size {}.'.format(final_dimension))

        o_shape = array.shape

        if array.ndim != 2:
            array = numpy.reshape(array, (-1, final_dimension))
        return array, o_shape

    def ecf_to_ortho(self, coords):
        """
        Gets the `(ortho_row, ortho_column)` coordinates in the ortho-rectified
        system for the provided physical coordinates in ECF `(X, Y, Z)` coordinates.

        Parameters
        ----------
        coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def ecf_to_pixel(self, coords):
        """
        Gets the `(pixel_row, pixel_column)` coordinates for the provided physical
        coordinates in ECF `(X, Y, Z)` coordinates.

        Parameters
        ----------
        coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def ll_to_ortho(self, ll_coords):
        """
        Gets the `(ortho_row, ortho_column)` coordinates in the ortho-rectified
        system for the provided physical coordinates in `(Lat, Lon)` coordinates.

        Note that there is inherent ambiguity when handling the missing elevation,
        and the effect is likely methodology dependent.

        Parameters
        ----------
        ll_coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def llh_to_ortho(self, llh_coords):
        """
        Gets the `(ortho_row, ortho_column)` coordinates in the ortho-rectified
        system for the providednphysical coordinates in `(Lat, Lon, HAE)`
        coordinates.

        Parameters
        ----------
        llh_coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def pixel_to_ortho(self, pixel_coords):
        """
        Gets the ortho-rectified indices for the point(s) in pixel coordinates.

        Parameters
        ----------
        pixel_coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def pixel_to_ecf(self, pixel_coords):
        """
        Gets the ECF coordinates for the point(s) in pixel coordinates.

        Parameters
        ----------
        pixel_coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def ortho_to_ecf(self, ortho_coords):
        """
        Get the ecf coordinates for the point(s) in ortho-rectified coordinates.

        Parameters
        ----------
        ortho_coords : numpy.ndarray
            Point(s) in the ortho-recitified coordinate system, of the form
            `(ortho_row, ortho_column)`.

        Returns
        -------
        numpy.ndarray
        """

        raise NotImplementedError

    def ortho_to_llh(self, ortho_coords):
        """
        Get the lat/lon/hae coordinates for the point(s) in ortho-rectified coordinates.

        Parameters
        ----------
        ortho_coords : numpy.ndarray
            Point(s) in the ortho-recitified coordinate system, of the form
            `(ortho_row, ortho_column)`.

        Returns
        -------
        numpy.ndarray
        """

        ecf = self.ortho_to_ecf(ortho_coords)
        return ecf_to_geodetic(ecf)

    def ortho_to_pixel(self, ortho_coords):
        """
        Get the pixel indices for the point(s) in ortho-rectified coordinates.

        Parameters
        ----------
        ortho_coords : numpy.ndarray
            Point(s) in the ortho-recitified coordinate system, of the form
            `(ortho_row, ortho_column)`.

        Returns
        -------
        numpy.ndarray
            The array of indices, of the same shape as `new_coords`, which indicate
            `(row, column)` pixel (fractional) indices.
        """

        raise NotImplementedError

    def get_pixel_array_bounds(self, coords):
        """
        Extract integer bounds of the input array, expected to have final dimension
        of size 2.

        Parameters
        ----------
        coords : numpy.ndarray

        Returns
        -------
        numpy.ndarray
            Of the form `(min_row, max_row, min_column, max_column)`.
        """

        coords, o_shape = self._reshape(coords, 2)
        return numpy.array(
            (numpy.ceil(numpy.nanmin(coords[:, 0], axis=0)),
             numpy.floor(numpy.nanmax(coords[:, 0], axis=0)),
             numpy.ceil(numpy.nanmin(coords[:, 1], axis=0)),
             numpy.floor(numpy.nanmax(coords[:, 1], axis=0))), dtype=numpy.int64)

class PGProjection(ProjectionHelper):
    """
    Class which helps perform the Planar Grid (i.e. Ground Plane) ortho-rectification
    for a sicd-type object using the SICD projection model directly.
    """

    __slots__ = (
        '_reference_point', '_reference_pixels', '_row_vector', '_col_vector', '_normal_vector', '_reference_hae')

    def __init__(self, sicd, reference_point=None, reference_pixels=None, normal_vector=None, row_vector=None,
                 col_vector=None, row_spacing=None, col_spacing=None,
                 default_pixel_method='GEOM_MEAN'):
        r"""

        Parameters
        ----------
        sicd : SICDType
            The sicd object
        reference_point : None|numpy.ndarray
            The reference point (origin) of the planar grid. If None, a default
            derived from the SICD will be used.
        reference_pixels : None|numpy.ndarray
            The projected pixel
        normal_vector : None|numpy.ndarray
            The unit normal vector of the plane.
        row_vector : None|numpy.ndarray
            The vector defining increasing column direction. If None, a default
            derived from the SICD will be used.
        col_vector : None|numpy.ndarray
            The vector defining increasing column direction. If None, a default
            derived from the SICD will be used.
        row_spacing : None|float
            The row pixel spacing.
        col_spacing : None|float
            The column pixel spacing.
        default_pixel_method : str
            Must be one of ('MAX', 'MIN', 'MEAN', 'GEOM_MEAN'). This determines
            the default behavior for row_spacing/col_spacing. The default value for
            row/column spacing will be the implied function applied to the range
            and azimuth ground resolution. Note that geometric mean is defined as
            :math:`\sqrt(x*x + y*y)`
        """

        self._reference_point = None
        self._reference_hae = None
        self._reference_pixels = None
        self._normal_vector = None
        self._row_vector = None
        self._col_vector = None
        ProjectionHelper.__init__(
            self, sicd, row_spacing=row_spacing, col_spacing=col_spacing, default_pixel_method=default_pixel_method)
        self.set_reference_point(reference_point=reference_point)
        self.set_reference_pixels(reference_pixels=reference_pixels)
        self.set_plane_frame(
            normal_vector=normal_vector, row_vector=row_vector, col_vector=col_vector)

    @property
    def reference_point(self):
        # type: () -> numpy.ndarray
        """
        numpy.ndarray: The grid reference point.
        """

        return self._reference_point

    @property
    def reference_pixels(self):
        # type: () -> numpy.ndarray
        """
        numpy.ndarray: The ortho-rectified pixel coordinates of the grid reference point.
        """

        return self._reference_pixels

    @property
    def normal_vector(self):
        # type: () -> numpy.ndarray
        """
        numpy.ndarray: The normal vector.
        """

        return self._normal_vector

    def set_reference_point(self, reference_point=None):
        """
        Sets the reference point, which must be provided in ECF coordinates.

        Parameters
        ----------
        reference_point : None|numpy.ndarray
            The reference point (origin) of the planar grid. If None, then the
            `sicd.GeoData.SCP.ECF` will be used.

        Returns
        -------
        None
        """

        if reference_point is None:
            if self.sicd.RadarCollection.Area is None:
                reference_point = self.sicd.GeoData.SCP.ECF.get_array()
            else:
                reference_point = self.sicd.RadarCollection.Area.Plane.RefPt.ECF.get_array()

        if not (isinstance(reference_point, numpy.ndarray) and reference_point.ndim == 1
                and reference_point.size == 3):
            raise ValueError('reference_point must be a vector of size 3.')
        self._reference_point = reference_point
        # set the reference hae
        ref_llh = ecf_to_geodetic(reference_point)
        self._reference_hae = ref_llh[2]

    def set_reference_pixels(self, reference_pixels=None):
        """
        Sets the reference point, which must be provided in ECF coordinates.

        Parameters
        ----------
        reference_pixels : None|numpy.ndarray
            The ortho-rectified pixel coordinates for the reference point (origin) of the planar grid.
            If None, then the (0, 0) will be used.

        Returns
        -------
        None
        """

        if reference_pixels is None:
            if self.sicd.RadarCollection.Area is not None:
                reference_pixels = numpy.array([
                    self.sicd.RadarCollection.Area.Plane.RefPt.Line,
                    self.sicd.RadarCollection.Area.Plane.RefPt.Sample],
                    dtype='float64')
            else:
                reference_pixels = numpy.zeros((2, ), dtype='float64')

        if not (isinstance(reference_pixels, numpy.ndarray) and reference_pixels.ndim == 1
                and reference_pixels.size == 2):
            raise ValueError('reference_pixels must be a vector of size 2.')
        self._reference_pixels = reference_pixels

    @property
    def row_vector(self):
        """
        numpy.ndarray: The grid increasing row direction (ECF) unit vector.
        """

        return self._row_vector

    @property
    def col_vector(self):
        """
        numpy.ndarray: The grid increasing column direction (ECF) unit vector.
        """

        return self._col_vector

    @property
    def reference_hae(self):
        """
        float: The height above the ellipsoid of the reference point.
        """

        return self._reference_hae

    def set_plane_frame(self, normal_vector=None, row_vector=None, col_vector=None):
        """
        Set the plane unit normal, and the row and column vectors, in ECF coordinates.
        Note that the perpendicular component of col_vector with respect to the
        row_vector will be used.

        If `normal_vector`, `row_vector`, and `col_vector` are all `None`, then
        the normal to the Earth tangent plane at the reference point is used for
        `normal_vector`. The `row_vector` will be defined as the perpendicular
        component of `sicd.Grid.Row.UVectECF` to `normal_vector`. The `colummn_vector`
        will be defined as the component of `sicd.Grid.Col.UVectECF` perpendicular
        to both `normal_vector` and `row_vector`.

        If only `normal_vector` is supplied, then the `row_vector` and `column_vector`
        will be defined similarly as the perpendicular components of
        `sicd.Grid.Row.UVectECF` and `sicd.Grid.Col.UVectECF`.

        Otherwise, all vectors supplied will be normalized, but are required to be
        mutually perpendicular. If only two vectors are supplied, then the third
        will be determined.

        Parameters
        ----------
        normal_vector : None|numpy.ndarray
            The vector defining the outward unit normal in ECF coordinates.
        row_vector : None|numpy.ndarray
            The vector defining increasing column direction.
        col_vector : None|numpy.ndarray
            The vector defining increasing column direction.

        Returns
        -------
        None
        """

        def normalize(vec, name, perp=None):
            if not isinstance(vec, numpy.ndarray):
                vec = numpy.array(vec, dtype=numpy.float64)
            if not (isinstance(vec, numpy.ndarray) and vec.ndim == 1 and vec.size == 3):
                raise ValueError('{} vector must be a numpy.ndarray of dimension 1 and size 3.'.format(name))
            vec = numpy.copy(vec)
            if perp is None:
                pass
            elif isinstance(perp, numpy.ndarray):
                vec = vec - perp*(perp.dot(vec))
            else:
                for entry in perp:
                    vec = vec - entry*(entry.dot(vec))

            norm = numpy.linalg.norm(vec)
            if norm == 0:
                raise ValueError('{} vector cannot be the zero vector.'.format(name))
            elif norm != 1:
                vec = vec/norm  # avoid modifying row_vector def exterior to this class
            return vec

        def check_perp(vec1, vec2, name1, name2, tolerance=1e-6):
            if abs(vec1.dot(vec2)) > tolerance:
                raise ValueError('{} vector and {} vector are not perpendicular'.format(name1, name2))

        if self._reference_point is None:
            raise ValueError('This requires that reference point is previously set.')

        if normal_vector is None and row_vector is None and col_vector is None:
            if self.sicd.RadarCollection.Area is None:
                self._normal_vector = wgs_84_norm(self.reference_point)
                self._row_vector = normalize(
                    self.sicd.Grid.Row.UVectECF.get_array(), 'row', perp=self.normal_vector)
                self._col_vector = normalize(
                    self.sicd.Grid.Col.UVectECF.get_array(), 'column', perp=(self.normal_vector, self.row_vector))
            else:
                self._row_vector = self.sicd.RadarCollection.Area.Plane.XDir.UVectECF.get_array()
                self._col_vector = normalize(
                    self.sicd.RadarCollection.Area.Plane.YDir.UVectECF.get_array(), 'col', perp=self._row_vector)
                self._normal_vector = numpy.cross(self._row_vector, self._col_vector)
        elif normal_vector is not None and row_vector is None and col_vector is None:
            self._normal_vector = normalize(normal_vector, 'normal')
            self._row_vector = normalize(
                self.sicd.Grid.Row.UVectECF.get_array(), 'row', perp=self.normal_vector)
            self._col_vector = normalize(
                self.sicd.Grid.Col.UVectECF.get_array(), 'column', perp=(self.normal_vector, self.row_vector))
        elif normal_vector is None:
            if row_vector is None or col_vector is None:
                raise ValueError('normal_vector is not defined, so both row_vector and col_vector must be.')
            row_vector = normalize(row_vector, 'row')
            col_vector = normalize(col_vector, 'col')
            check_perp(row_vector, col_vector, 'row', 'col')
            self._row_vector = row_vector
            self._col_vector = col_vector
            self._normal_vector = numpy.cross(row_vector, col_vector)
        elif col_vector is None:
            if row_vector is None:
                raise ValueError('col_vector is not defined, so both normal_vector and row_vector must be.')
            normal_vector = normalize(normal_vector, 'normal')
            row_vector = normalize(row_vector, 'row')
            check_perp(normal_vector, row_vector, 'normal', 'row')
            self._normal_vector = normal_vector
            self._row_vector = row_vector
            self._col_vector = numpy.cross(self.normal_vector, self.row_vector)
        elif row_vector is None:
            normal_vector = normalize(normal_vector, 'normal')
            col_vector = normalize(col_vector, 'col')
            check_perp(normal_vector, col_vector, 'normal', 'col')
            self._normal_vector = normal_vector
            self._col_vector = col_vector
            self._row_vector = numpy.cross(self.col_vector, self.normal_vector)
        else:
            normal_vector = normalize(normal_vector, 'normal')
            row_vector = normalize(row_vector, 'row')
            col_vector = normalize(col_vector, 'col')
            check_perp(normal_vector, row_vector, 'normal', 'row')
            check_perp(normal_vector, col_vector, 'normal', 'col')
            check_perp(row_vector, col_vector, 'row', 'col')
            self._normal_vector = normal_vector
            self._row_vector = row_vector
            self._col_vector = col_vector
        # check for outward unit norm
        if numpy.dot(self.normal_vector, self.reference_point) < 0:
            logger.warning(
                'The normal vector appears to be outward pointing, so reversing.')
            self._normal_vector *= -1

    def plane_ecf_to_ortho(self, coords):
        """
        Converts ECF coordinates **known to be in the ground plane** to ortho grid coordinates.

        Parameters
        ----------
        coords : numpy.ndarray

        Returns
        -------
        numpy.ndarray
        """

        coords, o_shape = self._reshape(coords, 3)
        diff = coords - self.reference_point
        if len(o_shape) == 1:
            out = numpy.zeros((2, ), dtype=numpy.float64)
            out[0] = self._reference_pixels[0] + numpy.sum(diff*self.row_vector)/self.row_spacing
            out[1] = self._reference_pixels[1] + numpy.sum(diff*self.col_vector)/self.col_spacing
        else:
            out = numpy.zeros((coords.shape[0], 2), dtype=numpy.float64)
            out[:, 0] = self._reference_pixels[0] + numpy.sum(diff*self.row_vector, axis=1)/self.row_spacing
            out[:, 1] = self._reference_pixels[1] + numpy.sum(diff*self.col_vector, axis=1)/self.col_spacing
            out = numpy.reshape(out, o_shape[:-1] + (2, ))
        return out

    def ecf_to_ortho(self, coords):
        return self.pixel_to_ortho(self.ecf_to_pixel(coords))

    def ecf_to_pixel(self, coords):
        pixel, _, _ = self.sicd.project_ground_to_image(coords, tolerance=1e-6, max_iterations=40)
        return pixel

    def ll_to_ortho(self, ll_coords):
        """
        Gets the `(ortho_row, ortho_column)` coordinates in the ortho-rectified
        system for the provided physical coordinates in `(Lat, Lon)` coordinates.
        In this case, the missing altitude will be set to `reference_hae`, which
        is imperfect.

        Parameters
        ----------
        ll_coords : numpy.ndarray|list|tuple

        Returns
        -------
        numpy.ndarray
        """

        ll_coords, o_shape = self._reshape(ll_coords, 2)
        llh_temp = numpy.zeros((ll_coords.shape[0], 3), dtype=numpy.float64)
        llh_temp[:, :2] = ll_coords
        llh_temp[:, 2] = self.reference_hae
        llh_temp = numpy.reshape(llh_temp, o_shape[:-1] + (3, ))
        return self.llh_to_ortho(llh_temp)

    def llh_to_ortho(self, llh_coords):
        llh_coords, o_shape = self._reshape(llh_coords, 3)
        ground = geodetic_to_ecf(llh_coords)
        return self.ecf_to_ortho(numpy.reshape(ground, o_shape))

    def ortho_to_ecf(self, ortho_coords):
        ortho_coords, o_shape = self._reshape(ortho_coords, 2)
        xs = (ortho_coords[:, 0] - self._reference_pixels[0])*self.row_spacing
        ys = (ortho_coords[:, 1] - self._reference_pixels[1])*self.col_spacing
        if xs.ndim == 0:
            coords = self.reference_point + xs*self.row_vector + ys*self.col_vector
        else:
            coords = self.reference_point + numpy.outer(xs, self.row_vector) + \
                     numpy.outer(ys, self.col_vector)
        return numpy.reshape(coords, o_shape[:-1] + (3, ))

    def ortho_to_pixel(self, ortho_coords):
        ortho_coords, o_shape = self._reshape(ortho_coords, 2)
        pixel, _, _ = self.sicd.project_ground_to_image(self.ortho_to_ecf(ortho_coords), tolerance=1e-3, max_iterations=25)
        return numpy.reshape(pixel, o_shape)

    def pixel_to_ortho(self, pixel_coords):
        return self.plane_ecf_to_ortho(self.pixel_to_ecf(pixel_coords))

    def pixel_to_ecf(self, pixel_coords):
        return self.sicd.project_image_to_ground(
            pixel_coords, projection_type='PLANE',
            gref=self.reference_point, ugpn=self.normal_vector)

class PGRatPolyProjection(PGProjection):

    __slots__ = (
        '_reference_point', '_reference_pixels', '_row_vector', '_col_vector',
        '_normal_vector', '_reference_hae',
        '_row_samples', '_col_samples', '_alt_samples', '_alt_span',
        '_ecf_to_pixel_func', '_pixel_to_ortho_func', '_ortho_to_pixel_func')

    def __init__(self, sicd, reference_point=None, reference_pixels=None, normal_vector=None, row_vector=None,
                 col_vector=None, row_spacing=None, col_spacing=None,
                 default_pixel_method='GEOM_MEAN',
                 row_samples=51, col_samples=51, alt_samples=11, alt_span=250):
        r"""

        Parameters
        ----------
        sicd : SICDType
            The sicd object
        reference_point : None|numpy.ndarray
            The reference point (origin) of the planar grid. If None, a default
            derived from the SICD will be used.
        reference_pixels : None|numpy.ndarray
            The projected pixel
        normal_vector : None|numpy.ndarray
            The unit normal vector of the plane.
        row_vector : None|numpy.ndarray
            The vector defining increasing column direction. If None, a default
            derived from the SICD will be used.
        col_vector : None|numpy.ndarray
            The vector defining increasing column direction. If None, a default
            derived from the SICD will be used.
        row_spacing : None|float
            The row pixel spacing.
        col_spacing : None|float
            The column pixel spacing.
        default_pixel_method : str
            Must be one of ('MAX', 'MIN', 'MEAN', 'GEOM_MEAN'). This determines
            the default behavior for row_spacing/col_spacing. The default value for
            row/column spacing will be the implied function applied to the range
            and azimuth ground resolution. Note that geometric mean is defined as
            :math:`\sqrt(x*x + y*y)`
        row_samples : int
            How many row samples to use in fitting
        col_samples : int
            How many column samples to use in fitting
        alt_samples : int
            How many altitude samples to use in fitting
        alt_span : int|float
            Fitting for reference point hae +/- alt_span information.
        """

        self._ecf_to_pixel_func = None
        self._pixel_to_ortho_func = None
        self._ortho_to_pixel_func = None
        self._row_samples = int(row_samples)
        self._col_samples = int(col_samples)
        self._alt_samples = int(alt_samples)
        self._alt_span = float(alt_span)

        PGProjection.__init__(
            self, sicd, reference_point=reference_point, reference_pixels=reference_pixels,
            normal_vector=normal_vector, row_vector=row_vector, col_vector=col_vector,
            row_spacing=row_spacing, col_spacing=col_spacing, default_pixel_method=default_pixel_method)
        self.perform_rational_poly_fitting()

    def _perform_ecf_func_fitting(self):
        num_rows = self.sicd.ImageData.NumRows
        num_cols = self.sicd.ImageData.NumCols

        row_array = numpy.linspace(0, num_rows-1, self._row_samples)
        col_array = numpy.linspace(0, num_cols-1, self._col_samples)
        hae_array = self.reference_hae + numpy.linspace(
            -self._alt_span, self._alt_span, self._alt_samples)

        row_col_grid = numpy.empty((row_array.size, col_array.size, 2), dtype='float64')
        row_col_grid[:, :, 1], row_col_grid[:, :, 0] = numpy.meshgrid(col_array, row_array)

        ECF_data = numpy.empty((row_array.size, col_array.size, hae_array.size, 3))
        for i, hae0 in enumerate(hae_array):
            ECF_data[:, :, i, :] = self.sicd.project_image_to_ground(row_col_grid, projection_type='HAE', hae0=hae0)

        if not numpy.all(numpy.isfinite(ECF_data)):
            raise SarpyRatPolyError(
                'NaN values are encountered when projecting across the image area,\n\t'
                'this SICD is not a good candidate for projection using rational polynomials')
        row_func = get_rational_poly_3d(
            ECF_data[:, :, :, 0].flatten(), ECF_data[:, :, :, 1].flatten(), ECF_data[:, :, :, 2].flatten(),
            numpy.stack([row_col_grid[:, :, 0] for _ in range(self._alt_samples)], axis=2).flatten(), order=3)
        col_func = get_rational_poly_3d(
            ECF_data[:, :, :, 0].flatten(), ECF_data[:, :, :, 1].flatten(), ECF_data[:, :, :, 2].flatten(),
            numpy.stack([row_col_grid[:, :, 1] for _ in range(self._alt_samples)], axis=2).flatten(), order=3)
        self._ecf_to_pixel_func = CombinedRationalPolynomial(row_func, col_func)

    def _perform_pixel_fitting(self):
        num_rows = self.sicd.ImageData.NumRows
        num_cols = self.sicd.ImageData.NumCols
        row_array = numpy.linspace(0, num_rows-1, 2*self._row_samples)
        col_array = numpy.linspace(0, num_cols-1, 2*self._col_samples)
        pixel_data = numpy.empty((row_array.size, col_array.size, 2), dtype='float64')
        pixel_data[:, :, 1], pixel_data[:, :, 0] = numpy.meshgrid(col_array, row_array)

        ortho_data = PGProjection.pixel_to_ortho(self, pixel_data)
        pix_to_orth_row = get_rational_poly_2d(
            pixel_data[:, :, 0].flatten(), pixel_data[:, :, 1].flatten(),
            ortho_data[:, :, 0], order=3)
        pix_to_orth_col = get_rational_poly_2d(
            pixel_data[:, :, 0].flatten(), pixel_data[:, :, 1].flatten(),
            ortho_data[:, :, 1], order=3)
        self._pixel_to_ortho_func = CombinedRationalPolynomial(pix_to_orth_row, pix_to_orth_col)

        orth_to_pix_row = get_rational_poly_2d(
            ortho_data[:, :, 0].flatten(), ortho_data[:, :, 1].flatten(),
            pixel_data[:, :, 0], order=3)
        orth_to_pix_col = get_rational_poly_2d(
            ortho_data[:, :, 0].flatten(), ortho_data[:, :, 1].flatten(),
            pixel_data[:, :, 1], order=3)
        self._ortho_to_pixel_func = CombinedRationalPolynomial(orth_to_pix_row, orth_to_pix_col)

    def perform_rational_poly_fitting(self):
        """
        Defined the rational polynomial functions via fitting.
        """

        self._perform_ecf_func_fitting()
        self._perform_pixel_fitting()

    def ecf_to_ortho(self, coords):
        return self.pixel_to_ortho(self.ecf_to_pixel(coords))

    def ecf_to_pixel(self, coords):
        return self._ecf_to_pixel_func(coords, combine=True)

    def ortho_to_pixel(self, ortho_coords):
        return self._ortho_to_pixel_func(ortho_coords, combine=True)

    def pixel_to_ortho(self, pixel_coords):
        return self._pixel_to_ortho_func(pixel_coords, combine=True)

max_stars_repo_name: spacefan/sarpy

max_issues_repo_name: spacefan/sarpy

max_forks_repo_name: spacefan/sarpy


text:
import pandas as pd
import numpy as np

# number of checkpoints
num_required = 10

def get_non_dom_ind(vals_arr, verbose=True):
    res_arr = []
    for i in range(0, len(vals_arr)):
        point1 = vals_arr[i]
        is_dominated = False
        for j in range(0, len(vals_arr)):
            if i == j:
                continue
            point2 = vals_arr[j]
            if ((point2[0] < point1[0]) and (point2[1] <= point1[1])) or (
                    (point2[0] <= point1[0]) and (point2[1] < point1[1])):
                if verbose:
                    print('{} is dominated by {}'.format(i, j))
                is_dominated = True
                break
            pass
        if not is_dominated:
            res_arr.append(i)
        pass
    return res_arr

def get_best_n_points_no_curve(n, arr_input):
    arr_list = arr_input.tolist()
    arr_list.sort(reverse=True, key=lambda x: x[0])
    arr_list = np.array(arr_list)

    # 2. find arc length
    arc_len_arr = []
    for pos in range(0, len(arr_input) - 1):
        p1 = np.array([arr_list[pos][0], arr_list[pos][1]])
        p2 = np.array([arr_list[pos + 1][0], arr_list[pos + 1][1]])
        arc_len_arr.append(np.linalg.norm(p2 - p1))
    arc_len_arr = np.array(arc_len_arr)
    # distance delta
    d = sum(arc_len_arr) / (n - 1)
    # cumul_sum of art length
    arc_len_arr_cum = np.cumsum(arc_len_arr)

    # 3. choose ref. points
    # positions of reference points
    points_pos = [0]
    j = 1
    for i in range(0, len(arc_len_arr_cum)):
        if arc_len_arr_cum[i] >= j * d:
            points_pos.append(i + 1)
            j += 1
            if j == n - 1:
                break
        pass
    points_pos.append(len(arr_list) - 1)

    chosen_points = []
    for ref_point_pos in points_pos:
        ref_point = arr_list[ref_point_pos]
        dist = np.linalg.norm((arr_input - ref_point), axis=1)
        pos = np.argmin(dist)
        chosen_points.append(pos)
        pass

    return chosen_points

def aggregate_results(results):
    df_arr = []
    for el in results:
        df = pd.read_csv(el)
        df.name = el
        df_arr.append(df)

    # get all test
    all_tests = []
    for df in df_arr:
        all_tests.extend(df['TestName'].unique())
        pass
    # remove duplicates
    all_tests = list(set(all_tests))

    # test = all_tests[0]

    final_dfs = []
    for test in all_tests:
        print('Processing {}'.format(test))
        # generate sub_dfs that will contain only the results for the required test
        test_df_arr = []
        for df in df_arr:
            test_df = df.loc[df['TestName'] == test].copy()
            test_df.name = df.name
            test_df['Geo'] = test_df['CountryName'].astype(str) + '__' + test_df['RegionName'].astype(str)
            test_df_arr.append(test_df)
            pass

        # get all geo-s in this test
        all_geos = set()
        for test_df in test_df_arr:
            geos = test_df['Geo'].unique()
            all_geos.update(set(geos))
            pass
        all_geos = list(all_geos)

        # for every geo
        for curr_geo in all_geos:
            # put all prescriptors for the current geo in the same df
            prescr_df = []
            for test_df in test_df_arr:
                df = test_df.loc[test_df['Geo'] == curr_geo].copy()
                df['source'] = test_df.name
                prescr_df.append(df)
                pass
            # put them all together and remove duplicates
            prescr_df = pd.concat(prescr_df).drop_duplicates(['PredictedDailyNewCases', 'Stringency'])

            # choose non-dominated
            arr = prescr_df[['PredictedDailyNewCases', 'Stringency']].values
            non_dom_idx = get_non_dom_ind(arr, verbose=False)

            # check here if we have at least 10
            if len(non_dom_idx) < num_required:
                # print('less than 10: {}'.format(curr_geo))
                # just copy the first point the required number of times
                chosen_points_pos = non_dom_idx
                num_found = len(chosen_points_pos)
                for i in range(0, num_required - num_found):
                    chosen_points_pos.append(chosen_points_pos[0])
                    pass
                chosen_non_dom = chosen_points_pos
                pass
            elif len(non_dom_idx) == num_required:
                # don't do anything
                chosen_points_pos = non_dom_idx
                chosen_non_dom = chosen_points_pos
                pass
            else:
                # choose points
                chosen_points_pos = get_best_n_points_no_curve(num_required, arr[non_dom_idx])
                chosen_non_dom = np.array(non_dom_idx)[chosen_points_pos]
                pass
            # create a df with chosen tests for this geo
            geo_df = prescr_df.iloc[chosen_non_dom].copy()
            geo_df.drop('Unnamed: 0', axis=1, inplace=True)
            tmp = geo_df.copy()
            while len(geo_df) < num_required:
                geo_df = pd.concat([geo_df, tmp])
            geo_df = geo_df[0: num_required]
            geo_df['source-PrescriptionIndex'] = geo_df['PrescriptionIndex']
            geo_df['PrescriptionIndex'] = [i for i in range(0, num_required)]
            final_dfs.append(geo_df)
            pass
        pass
    res_df = pd.concat(final_dfs)
    res_df.reset_index(drop=True, inplace=True)
    return res_df

# path to files with results
results = [
        'Heuristic_pct_5.csv',
        'neat_evaluate_2D_big_test.csv',
    ]
final_df = aggregate_results(results)

final_df.to_csv('aggregation.csv')

max_stars_repo_name: bradyneal/covid-xprize-comp

max_issues_repo_name: bradyneal/covid-xprize-comp

max_forks_repo_name: bradyneal/covid-xprize-comp


text:
from collections import OrderedDict
import os
import time
import numpy as np
from tqdm import tqdm

from torch.autograd import Variable
import torch.optim as optim
import torch.nn as nn
import torch.utils.data
import torchvision.utils as vutils
from models.basemodel import *

from models.networks import NetG, NetD, weights_init
from models.evaluate import evaluate
from utils.visualizer import Visualizer
from utils.loss import l2_loss

import random

class BaseModel():
    """ Base Model for ganomaly
    """

    def __init__(self, opt):
        ##
        # Seed for deterministic behavior
        self.seed(opt.manualseed)

        # Initalize variables.
        self.opt = opt
        self.visualizer = Visualizer(opt)
        self.trn_dir = os.path.join(self.opt.outf, self.opt.name, 'train')
        self.tst_dir = os.path.join(self.opt.outf, self.opt.name, 'test')
        self.device = torch.device("cuda:0" if self.opt.device != 'cpu' else "cpu")
        self.dataloader_setup()

    def dataloader_setup(self):
        raise NotImplementedError("dataloader_setup not implemented")
    ##
    def set_input(self, input: torch.Tensor):
        """ Set input and ground truth

        Args:
            input (FloatTensor): Input data for batch i.
        """
        with torch.no_grad():
            self.input.resize_(input[0].size()).copy_(input[0])
            self.gt.resize_(input[1].size()).copy_(input[1])
            self.label.resize_(input[1].size())

            # Copy the first batch as the fixed input.
            if self.total_steps == self.opt.batchsize:
                self.fixed_input.resize_(input[0].size()).copy_(input[0])

    ##
    def seed(self, seed_value):
        """ Seed

        Arguments:
            seed_value {int} -- [description]
        """
        # Check if seed is default value
        if seed_value == -1:
            return

        # Otherwise seed all functionality
        import random
        random.seed(seed_value)
        torch.manual_seed(seed_value)
        torch.cuda.manual_seed_all(seed_value)
        np.random.seed(seed_value)
        torch.backends.cudnn.deterministic = True

    ##
    def get_errors(self):
        """ Get netD and netG errors.

        Returns:
            [OrderedDict]: Dictionary containing errors.
        """

        errors = OrderedDict([
            ('err_d', self.err_d.item()),
            ('err_g', self.err_g.item()),
            ('err_g_adv', self.err_g_adv.item()),
            ('err_g_con', self.err_g_con.item()),
            ('err_g_enc', self.err_g_enc.item())])

        return errors

    ##
    def get_current_images(self):
        """ Returns current images.

        Returns:
            [reals, fakes, fixed]
        """

        reals = self.input.data
        fakes = self.fake.data
        fixed = self.netgs[0](self.fixed_input)[0].data

        return reals, fakes, fixed

    ##
    def save_weights(self, epoch):
        """Save netG and netD weights for the current epoch.

        Args:
            epoch ([int]): Current epoch number.
        """

        weight_dir = os.path.join(self.opt.outf, self.opt.name, 'train', 'weights')
        if not os.path.exists(weight_dir): os.makedirs(weight_dir)
        for idx_g in range(self.opt.n_G):
            torch.save({'epoch': epoch + 1, 'state_dict': self.netgs[idx_g].state_dict()},
                       '%s/netG_%d.pth' % (weight_dir, idx_g))
        for idx_d in range(self.opt.n_D):
            torch.save({'epoch': epoch + 1, 'state_dict': self.netds[idx_d].state_dict()},
                       '%s/netD_%d.pth' % (weight_dir, idx_d))

    ##
    def train_one_epoch(self):
        """ Train the model for one epoch.
        """

        # self.netg.train()

        epoch_iter = 0
        for data in tqdm(self.dataloader.train, leave=False, total=len(self.dataloader.train)):
            self.total_steps += self.opt.batchsize
            epoch_iter += self.opt.batchsize

            idx_d = random.randint(0, self.opt.n_D - 1)
            idx_g = random.randint(0, self.opt.n_G - 1)

            self.netgs[idx_g].train()
            self.netds[idx_d].train()
            self.set_input(data)
            # self.optimize()
            self.optimize_params(idx_d, idx_g)

            if self.total_steps % self.opt.print_freq == 0:
                errors = self.get_errors()
                if self.opt.display:
                    counter_ratio = float(epoch_iter) / len(self.dataloader.train.dataset)
                    self.visualizer.plot_current_errors(self.epoch, counter_ratio, errors)

        print(">> Training model %s. Epoch %d/%d" % (self.name, self.epoch + 1, self.opt.niter))

    ##
    def train(self):
        """ Train the model
        """

        self.total_steps = 0
        best_auc = 0

        # Train for niter epochs.
        print(">> Training model %s." % self.name)
        for self.epoch in range(self.opt.iter, self.opt.niter):
            # Train for one epoch
            self.train_one_epoch()
            res = self.test()
            if res['AUC'] > best_auc:
                best_auc = res['AUC']
                self.save_weights(self.epoch)
            self.visualizer.print_current_performance(res, best_auc)
        print(">> Training model %s.[Done]" % self.name)

    ##
    def test(self):
        """ Test GANomaly model.

        Args:
            dataloader ([type]): Dataloader for the test set

        Raises:
            IOError: Model weights not found.
        """
        with torch.no_grad():
            # Load the weights of netg and netd.
            if self.opt.load_weights:
                path = "./output/{}/{}/train/weights/netG.pth".format(self.name.lower(), self.opt.dataset)
                pretrained_dict = torch.load(path)['state_dict']

                try:
                    self.netg.load_state_dict(pretrained_dict)
                except IOError:
                    raise IOError("netG weights not found")
                print('   Loaded weights.')

            self.opt.phase = 'test'

            # Create big error tensor for the test set.
            self.an_scores = torch.zeros(size=(len(self.dataloader.valid.dataset),), dtype=torch.float32,
                                         device=self.device)
            self.gt_labels = torch.zeros(size=(len(self.dataloader.valid.dataset),), dtype=torch.long,
                                         device=self.device)
            self.latent_i = torch.zeros(size=(len(self.dataloader.valid.dataset), self.opt.nz), dtype=torch.float32,
                                        device=self.device)
            self.latent_o = torch.zeros(size=(len(self.dataloader.valid.dataset), self.opt.nz), dtype=torch.float32,
                                        device=self.device)

            self.times = []
            self.total_steps = 0
            epoch_iter = 0
            for i, data in enumerate(self.dataloader.valid, 0):
                self.total_steps += self.opt.batchsize
                epoch_iter += self.opt.batchsize
                time_i = time.time()
                self.set_input(data)
                latent_is = []
                latent_os = []
                for idx_g in range(self.opt.n_G):
                    # self.fake, latent_i, latent_o = self.netgs[idx_g](self.input)
                    _, latent_i, latent_o = self.netgs[idx_g](self.input)
                    latent_is.append(latent_i)
                    latent_os.append(latent_o)

                latent_i = torch.mean(torch.stack(latent_is), dim=0)
                latent_o = torch.mean(torch.stack(latent_os), dim=0)

                error = torch.mean(torch.pow((latent_i - latent_o), 2), dim=1)
                time_o = time.time()

                self.an_scores[i * self.opt.batchsize: i * self.opt.batchsize + error.size(0)] = error.reshape(
                    error.size(0))
                self.gt_labels[i * self.opt.batchsize: i * self.opt.batchsize + error.size(0)] = self.gt.reshape(
                    error.size(0))
                self.latent_i[i * self.opt.batchsize: i * self.opt.batchsize + error.size(0), :] = latent_i.reshape(
                    error.size(0), self.opt.nz)
                self.latent_o[i * self.opt.batchsize: i * self.opt.batchsize + error.size(0), :] = latent_o.reshape(
                    error.size(0), self.opt.nz)

                self.times.append(time_o - time_i)

            # Measure inference time.
            self.times = np.array(self.times)
            self.times = np.mean(self.times[:100] * 1000)

            # Scale error vector between [0, 1]
            self.an_scores = (self.an_scores - torch.min(self.an_scores)) / (
                        torch.max(self.an_scores) - torch.min(self.an_scores))
            # auc, eer = roc(self.gt_labels, self.an_scores)
            auc = evaluate(self.gt_labels, self.an_scores, metric=self.opt.metric)
            performance = OrderedDict([('Avg Run Time (ms/batch)', self.times), ('AUC', auc)])

            if self.opt.display_id > 0 and self.opt.phase == 'test':
                counter_ratio = float(epoch_iter) / len(self.dataloader.valid.dataset)
                self.visualizer.plot_performance(self.epoch, counter_ratio, performance)
            return performance

##
class Ganomaly(BaseModel):
    """GANomaly Class
    """

    @property
    def name(self):
        return 'Ganomaly'

    def __init__(self, opt):
        super(Ganomaly, self).__init__(opt)

        # -- Misc attributes
        self.epoch = 0
        self.times = []
        self.total_steps = 0

        ##
        # Create and initialize networks.
        self.netgs = [NetG(self.opt).to(self.device) for _ in range(self.opt.n_G)]
        self.netds = [NetD(self.opt).to(self.device) for _ in range(self.opt.n_D)]
        for i in range(self.opt.n_G):
            self.netgs[i].apply(weights_init)
        for i in range(self.opt.n_D):
            self.netds[i].apply(weights_init)

        self.l_adv = l2_loss
        self.l_con = nn.L1Loss()
        self.l_enc = l2_loss
        self.l_bce = nn.BCELoss()

        ##
        # Initialize input tensors.
        self.input = torch.empty(size=(self.opt.batchsize, 3, self.opt.isize, self.opt.isize), dtype=torch.float32,
                                 device=self.device)
        self.label = torch.empty(size=(self.opt.batchsize,), dtype=torch.float32, device=self.device)
        self.gt = torch.empty(size=(opt.batchsize,), dtype=torch.long, device=self.device)
        self.fixed_input = torch.empty(size=(self.opt.batchsize, 3, self.opt.isize, self.opt.isize),
                                       dtype=torch.float32, device=self.device)
        self.real_label = torch.ones(size=(self.opt.batchsize,), dtype=torch.float32, device=self.device)
        self.fake_label = torch.zeros(size=(self.opt.batchsize,), dtype=torch.float32, device=self.device)
        ##
        # Setup optimizer
        if self.opt.isTrain:
            for i in range(self.opt.n_G):
                self.netgs[i].train()
            for i in range(self.opt.n_D):
                self.netds[i].train()
            self.optimizer_d = optim.Adam(sum([list(self.netds[i].parameters()) for i in range(self.opt.n_D)], []),
                                          lr=self.opt.lr, betas=(self.opt.beta1, 0.999))
            self.optimizer_g = optim.Adam(sum([list(self.netgs[i].parameters()) for i in range(self.opt.n_G)], []),
                                          lr=self.opt.lr, betas=(self.opt.beta1, 0.999))
    def dataloader_setup(self):
        self.dataloader = load_data(self.opt)

    ##
    def forward_g(self, idx_g):
        """ Forward propagate through netG
        """
        self.fake, self.latent_i, self.latent_o = self.netgs[idx_g](self.input)

    ##
    def forward_d(self, idx_d):
        """ Forward propagate through netD
        """
        self.pred_real, self.feat_real = self.netds[idx_d](self.input)
        self.pred_fake, self.feat_fake = self.netds[idx_d](self.fake.detach())

    ##
    def backward_g(self, idx_d):
        """ Backpropagate through netG
        """
        self.err_g_adv = self.l_adv(self.netds[idx_d](self.input)[1], self.netds[idx_d](self.fake)[1])
        self.err_g_con = self.l_con(self.fake, self.input)
        self.err_g_enc = self.l_enc(self.latent_o, self.latent_i)
        self.err_g = self.err_g_adv * self.opt.w_adv + \
                     self.err_g_con * self.opt.w_con + \
                     self.err_g_enc * self.opt.w_enc
        self.err_g.backward(retain_graph=True)

    ##
    def backward_d(self):
        """ Backpropagate through netD
        """
        # Real - Fake Loss
        self.err_d_real = self.l_bce(self.pred_real, self.real_label)
        self.err_d_fake = self.l_bce(self.pred_fake, self.fake_label)

        # NetD Loss & Backward-Pass
        self.err_d = (self.err_d_real + self.err_d_fake) * 0.5
        self.err_d.backward()

    ##
    def reinit_d(self, idx_d):
        """ Re-initialize the weights of netD
        """
        self.netds[idx_d].apply(weights_init)
        print('Reloading net d')

    def optimize_params(self, idx_d, idx_g):
        """ Forwardpass, Loss Computation and Backwardpass.
        """
        # Forward-pass
        self.forward_g(idx_g)
        self.forward_d(idx_d)

        # Backward-pass
        # netg
        self.optimizer_g.zero_grad()
        self.backward_g(idx_d)
        self.optimizer_g.step()

        # netd
        self.optimizer_d.zero_grad()
        self.backward_d()
        self.optimizer_d.step()
        # if self.err_d.item() < 1e-5: self.reinit_d(idx_d)


text:
# coding=utf-8
# Copyright 2020 HuggingFace Inc.
#
# 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.

import doctest
import glob
import importlib
import inspect
import os
from contextlib import contextmanager
from functools import wraps
from unittest.mock import patch

import numpy as np
from absl.testing import parameterized

import datasets
from datasets import load_metric

from .utils import for_all_test_methods, local, slow

REQUIRE_FAIRSEQ = {"comet"}
_has_fairseq = importlib.util.find_spec("fairseq") is not None

def skip_if_dataset_requires_fairseq(test_case):
    @wraps(test_case)
    def wrapper(self, metric_name):
        if not _has_fairseq and metric_name in REQUIRE_FAIRSEQ:
            self.skipTest('"test requires Fairseq"')
        else:
            test_case(self, metric_name)

    return wrapper

def get_local_metric_names():
    metrics = [metric_dir.split(os.sep)[-2] for metric_dir in glob.glob("./metrics/*/")]
    return [{"testcase_name": x, "metric_name": x} for x in metrics if x != "gleu"]  # gleu is unfinished

@parameterized.named_parameters(get_local_metric_names())
@for_all_test_methods(skip_if_dataset_requires_fairseq)
@local
class LocalMetricTest(parameterized.TestCase):
    INTENSIVE_CALLS_PATCHER = {}
    metric_name = None

    def test_load_metric(self, metric_name):
        doctest.ELLIPSIS_MARKER = "[...]"
        metric_module = importlib.import_module(datasets.load.prepare_module(os.path.join("metrics", metric_name))[0])
        metric = datasets.load.import_main_class(metric_module.__name__, dataset=False)
        # check parameters
        parameters = inspect.signature(metric._compute).parameters
        self.assertTrue("predictions" in parameters)
        self.assertTrue("references" in parameters)
        self.assertTrue(all([p.kind != p.VAR_KEYWORD for p in parameters.values()]))  # no **kwargs
        # run doctest
        with self.patch_intensive_calls(metric_name, metric_module.__name__):
            with self.use_local_metrics():
                results = doctest.testmod(metric_module, verbose=True, raise_on_error=True)
        self.assertEqual(results.failed, 0)
        self.assertGreater(results.attempted, 1)

    @slow
    def test_load_real_metric(self, metric_name):
        doctest.ELLIPSIS_MARKER = "[...]"
        metric_module = importlib.import_module(datasets.load.prepare_module(os.path.join("metrics", metric_name))[0])
        # run doctest
        with self.use_local_metrics():
            results = doctest.testmod(metric_module, verbose=True, raise_on_error=True)
        self.assertEqual(results.failed, 0)
        self.assertGreater(results.attempted, 1)

    @contextmanager
    def patch_intensive_calls(self, metric_name, module_name):
        if metric_name in self.INTENSIVE_CALLS_PATCHER:
            with self.INTENSIVE_CALLS_PATCHER[metric_name](module_name):
                yield
        else:
            yield

    @contextmanager
    def use_local_metrics(self):
        def load_local_metric(metric_name, *args, **kwargs):
            return load_metric(os.path.join("metrics", metric_name), *args, **kwargs)

        with patch("datasets.load_metric") as mock_load_metric:
            mock_load_metric.side_effect = load_local_metric
            yield

    @classmethod
    def register_intensive_calls_patcher(cls, metric_name):
        def wrapper(patcher):
            patcher = contextmanager(patcher)
            cls.INTENSIVE_CALLS_PATCHER[metric_name] = patcher
            return patcher

        return wrapper

# Metrics intensive calls patchers
# --------------------------------

@LocalMetricTest.register_intensive_calls_patcher("bleurt")
def patch_bleurt(module_name):
    import tensorflow.compat.v1 as tf
    from bleurt.score import Predictor

    tf.flags.DEFINE_string("sv", "", "")  # handle pytest cli flags

    class MockedPredictor(Predictor):
        def predict(self, input_dict):
            assert len(input_dict["input_ids"]) == 2
            return np.array([1.03, 1.04])

    # mock predict_fn which is supposed to do a forward pass with a bleurt model
    with patch("bleurt.score._create_predictor") as mock_create_predictor:
        mock_create_predictor.return_value = MockedPredictor()
        yield

@LocalMetricTest.register_intensive_calls_patcher("bertscore")
def patch_bertscore(module_name):
    import torch

    def bert_cos_score_idf(model, refs, *args, **kwargs):
        return torch.tensor([[1.0, 1.0, 1.0]] * len(refs))

    # mock get_model which is supposed to do download a bert model
    # mock bert_cos_score_idf which is supposed to do a forward pass with a bert model
    with patch("bert_score.scorer.get_model"), patch(
        "bert_score.scorer.bert_cos_score_idf"
    ) as mock_bert_cos_score_idf:
        mock_bert_cos_score_idf.side_effect = bert_cos_score_idf
        yield

@LocalMetricTest.register_intensive_calls_patcher("comet")
def patch_comet(module_name):
    def download_model(model):
        class Model:
            def predict(self, data, *args, **kwargs):
                assert len(data) == 2
                scores = [0.19, 0.92]
                data[0]["predicted_score"] = scores[0]
                data[1]["predicted_score"] = scores[1]
                return data, scores

        print("Download succeeded. Loading model...")
        return Model()

    # mock download_model which is supposed to do download a bert model
    with patch("comet.models.download_model") as mock_download_model:
        mock_download_model.side_effect = download_model
        yield

max_stars_repo_name: alexwdong/datasets

max_issues_repo_name: alexwdong/datasets

max_forks_repo_name: alexwdong/datasets


text:
import torch.nn as nn
import torch
import cv2
# x = torch.randn(60,1,3,3)
# # conv = torch.nn.Conv2d(1,64,3,padding=1)
# max = nn.MaxPool2d(2)
#
# res = max(x)
#
# print(res.shape)

ff = r"D:\datasets\dataaug\data\vitium\0\027_49.bmp 0"
a,v = ff.split(" ")
print(a)
print(v)
# import numpy as np
#
# a = torch.arange(0,6)
# b = a.view(2,3)
# c = b.unsqueeze(1).expand(2,2,3)
# print(c)

# import torch as t
# a = t.arange(0, 16).view(4, 4)
#
# print(a)
# # 选取对角线的元素
# index = t.LongTensor([[0, 1, 2, 3],[0, 1, 2, 3]])
# c = a.gather(0, index)
# print(c)

# t = torch.Tensor([[1, 2], [3,4]])
# print(t)
# g = torch.gather(t,1, torch.LongTensor([ [1, 0],[0,1]]))
# print(g)

a = torch.Tensor(range(0,30)).view(2,3,5)
print(a)

index = torch.LongTensor([[[0,1,2,0,2],
                           [0,0,0,0,0],
                           [1,1,1,1,1]],
                          [[1,2,2,2,2],
                           [0,0,0,0,0],
                           [2,2,2,2,2]]])
index2 = torch.LongTensor([[[0,1,1,0,1],
                            [0,1,1,1,1],
                            [1,1,1,1,1]],
                           [[1,0,0,0,0],
                            [0,0,0,0,0],
                            [1,1,0,0,0]],
                           [[1,0,0,0,0],
                            [0,0,0,0,0],
                            [1,1,0,0,0]]])
b = torch.gather(a,2,index)
print("dim=1:\n",b)

# import torch
# a = torch.randint(0, 30, (2, 3, 5))
# print(a)
#
# index2 = torch.LongTensor([[[0,1,1,0,1],
#                           [0,1,1,1,1],
#                           [1,1,1,1,1]],
#                         [[1,0,0,0,0],
#                          [0,0,0,0,0],
#                          [1,1,0,0,0]]])
# d = torch.gather(a, 0,index2)
# print(d)

img = cv2.imread(r"D:\project\prototypical_networks_master\data\omniglot\data\Anglo-Saxon_Futhorc\character04\0299_01.png")
img2 = cv2.imread(r"D:\datasets\dataaug\data\vitium\1\0149_32.bmp")

print(img.shape)
print(img2.shape)

print(int(4* 0.4))

max_stars_repo_name: loeweli/prototypical

max_issues_repo_name: loeweli/prototypical

max_forks_repo_name: loeweli/prototypical


text:
import numpy as np
import nrrd
import os
from os import makedirs
from os.path import join
from argparse import ArgumentParser
'''
Takes a folder full of .nrrd files and averages the lesion masks - outputs a "lesion prior".
Files that contain lesion masks should end in 'target.nrrd'. Lesion masks used in /usr/local/data/thomasc/unet_out/all_img
are linearly co-registered (MSLAQ).
TO DO - make a lesion prior with masks from non-linearly co-registered MSLAQ data (also
registered to patients)
'''

def main(args):
	in_dir = args.input
	out_dir = args.output
	makedirs(out_dir, exist_ok=True)

	#Only use training set for lesion prior
	subjects = np.load(in_dir)
	subjects = [i[4:] for i in subjects]
	nb_folds = 10
	fold_length = len(subjects) // nb_folds
	train_idx = (nb_folds - 2) * fold_length
	subjects = subjects[:train_idx]
	img_sum = 0
	img_count = 0

	for file in os.listdir(in_dir):
		if file.endswith("target.nrrd"):
			img, opts = nrrd.read(join(in_dir,file))
			img_sum = np.add(img_sum,img)
			img_count += 1
	print(img_count)
	les_prior = img_sum / img_count

	_OPTS = {'space': 'RAS', 'space directions': [(1, 0, 0), (0, 1, 0), (0, 0, 3)]}
	nrrd.write(join(out_dir, 'lesion_prior.nrrd'), les_prior, options=_OPTS)

def _parser():
	#Input path as /usr/local/data/tnair/thesis/data/mslaq/mslaq_subj_list.npy (blacktusk)
    usage = 'python les_prior.py -i input/path -o /output/path'
    parser = ArgumentParser()
    parser.add_argument('-i', '--input', help='Target Lesion Mask', required=True)
    parser.add_argument('-o', '--output', help='Lesion Prior Output', required=True)
    return parser

if __name__ == '__main__':
	main(_parser().parse_args())

max_forks_repo_name: thomaschristinck/lesion-detection-net


text:
"""A set of tools for use in integration tests."""
import os
import tempfile
import subprocess
import sys
from hashlib import sha1
import numpy as np
import tables
from nose.tools import assert_equal

CYCLUS_HAS_COIN = None

if sys.version_info[0] >= 3:
    str_types = (bytes, str)
else:
    str_types = (str, unicode)

def hasher(x):
    return int(sha1(x.encode()).hexdigest(), 16)

def idx(h):
    ind = [None] * 5
    for i in range(4, -1, -1):
        h, ind[i] = divmod(h, 2**32)
    return tuple(ind)

def sha1array(x): return np.array(idx(hasher(x)), np.uint32)

def table_exist(db, tables):
    """Checks if hdf5 database contains the specified tables.
    """
    return all([t in db.root for t in tables])

def find_ids(data, data_table, id_table):
    """Finds ids of the specified data located in the specified data_table,
    and extracts the corresponding id from the specified id_table.
    """
    ids = []
    for i, d in enumerate(data_table):
        if isinstance(d, np.ndarray) and isinstance(data, np.ndarray):
            if (d == data).all():
                ids.append(id_table[i])
        elif isinstance(d, np.ndarray) and not isinstance(data, np.ndarray):
            if (d == sha1array(data)).all():
                ids.append(id_table[i])
        elif d == data:
            ids.append(id_table[i])
    return ids

def exit_times(agent_id, exit_table):
    """Finds exit times of the specified agent from the exit table.
    """
    i = 0
    exit_times = []
    for index in exit_table["AgentId"]:
        if index == agent_id:
            exit_times.append(exit_table["ExitTime"][i])
        i += 1

    return exit_times

def run_cyclus(cyclus, cwd, in_path, out_path):
    """Runs cyclus with various inputs and creates output databases
    """
    holdsrtn = [1]  # needed because nose does not send() to test generator
    # make sure the output target directory exists
    cmd = [cyclus, "-o", out_path, "--input-file", in_path]
    check_cmd(cmd, cwd, holdsrtn)

def check_cmd(args, cwd, holdsrtn):
    """Runs a command in a subprocess and verifies that it executed properly.
    """
    if not isinstance(args, str_types):
        args = " ".join(args)
    print("TESTING: running command in {0}:\n\n{1}\n".format(cwd, args))
    env = dict(os.environ)
    env['_'] = subprocess.check_output(['which', 'cyclus'], cwd=cwd).strip()
    with tempfile.NamedTemporaryFile() as f:
        rtn = subprocess.call(args, shell=True, cwd=cwd, stdout=f, stderr=f,
                              env=env)
        if rtn != 0:
            f.seek(0)
            print("STDOUT + STDERR:\n\n" + f.read().decode())
    holdsrtn[0] = rtn
    assert_equal(rtn, 0)

def cyclus_has_coin():
    global CYCLUS_HAS_COIN
    if CYCLUS_HAS_COIN is not None:
        return CYCLUS_HAS_COIN
    s = subprocess.check_output(['cyclus', '--version'],
                                universal_newlines=True)
    s = s.strip().replace('Dependencies:', '')
    m = {k.strip(): v.strip() for k, v in [line.split()[:2] for line
                                           in s.splitlines() if line != '']}
    CYCLUS_HAS_COIN = m['Coin-Cbc'] != '-1'
    return CYCLUS_HAS_COIN


text:
import numpy as np
import tensorflow as tf
from tensorflow.keras import backend as K

def compute_accuracy(y_true, y_pred):
    '''Compute classification accuracy with a fixed threshold on distances.
    '''
    pred = y_pred.ravel() > 0.5
    return np.mean(pred == y_true)

def euclidean_distance(vects):
    x, y = vects
    sum_square = K.sum(K.square(x - y), axis=1, keepdims=True)
    return K.sqrt(K.maximum(sum_square, K.epsilon()))

def eucl_dist_output_shape(shapes):
    shape1, shape2 = shapes
    return (shape1[0], 1)

class ContrastiveLossWithMargin(tf.keras.losses.Loss):
    '''Contrastive loss from Hadsell-et-al.'06
    http://yann.lecun.com/exdb/publis/pdf/hadsell-chopra-lecun-06.pdf
    '''
    # initialize instance attributes
    def __init__(self, margin, name="contrastive_loss"):
        super().__init__(name=name)
        self.margin = margin

    # compute loss
    def call(self, y_true, y_pred):
        square_pred = K.square(y_pred)
        margin_square = K.square(K.maximum(self.margin - y_pred, 0))
        loss = K.mean(y_true * square_pred + (1 - y_true) * margin_square)
        return loss

max_stars_repo_name: tallamjr/dnn

max_issues_repo_name: tallamjr/dnn

max_forks_repo_name: tallamjr/dnn