text:
"""Manages the imaging context. This take a string and returns a dictionary containing:
 * Predict function
 * Invert function
 * image_iterator function
 * vis_iterator function

"""

import collections
import logging

import numpy

from data_models.memory_data_models import Image, Visibility
from data_models.parameters import get_parameter
from processing_library.image.operations import create_empty_image_like
from workflows.shared.imaging.imaging_shared import sum_invert_results, remove_sumwt, sum_predict_results
from wrappers.mpi.griddata.gridding import grid_weight_to_griddata, griddata_reweight, griddata_merge_weights
from wrappers.mpi.griddata.kernels import create_pswf_convolutionfunction
from wrappers.mpi.griddata.operations import create_griddata_from_image
from wrappers.mpi.image.deconvolution import deconvolve_cube, restore_cube
from wrappers.mpi.image.gather_scatter import image_gather_facets, \
    image_scatter_channels, image_gather_channels
from wrappers.mpi.image.operations import calculate_image_frequency_moments
from visibility import copy_visibility
from wrappers.mpi.imaging.weighting import taper_visibility_gaussian, taper_visibility_tukey

from mpi4py import MPI

log = logging.getLogger(__name__)

def predict_list_mpi_workflow(vis_list, model_imagelist, context, vis_slices=1, facets=1,
                                     gcfcf=None, comm=MPI.COMM_WORLD, **kwargs):
    """Predict, iterating over both the scattered vis_list and image

    The visibility and image are scattered, the visibility is predicted on each part, and then the
    parts are assembled. About data distribution: vis_list and model_imagelist
    live in rank 0; vis_slices, facets, context are replicated in all nodes.
    gcfcf if exists lives in rank 0; if not every mpi proc will create its own
    for the  corresponding subset of the image.

    :param vis_list:
    :param model_imagelist: Model used to determine image parameters
    :param vis_slices: Number of vis slices (w stack or timeslice)
    :param facets: Number of facets (per axis)
    :param context: Type of processing e.g. 2d, wstack, timeslice or facets
    :param gcfcg: tuple containing grid correction and convolution function
    :param comm: MPI communicator
    :param kwargs: Parameters for functions in components
    :return: List of vis_lists
   """
    rank = comm.Get_rank()
    size = comm.Get_size()
    log.info('%d: In predict_list_mpi_workflow: %d elements in vis_list' % (rank,len(vis_list)))
    # the assert only makes sense in proc 0 as for the others both lists are
    # empty
    assert len(vis_list) == len(model_imagelist), "Model must be the same length as the vis_list"

    # The use_serial_predict version paralelizes by freq (my opt version)
    assert get_parameter(kwargs, "use_serial_predict", True),"Only freq paralellization implemented"
    #if get_parameter(kwargs, "use_serial_predict", False):
    if get_parameter(kwargs, "use_serial_predict", True):
        from workflows.serial.imaging.imaging_serial import predict_list_serial_workflow

        image_results_list = list()
        # Distribute visibilities and model by freq
        sub_vis_list= numpy.array_split(vis_list, size)
        sub_vis_list=comm.scatter(sub_vis_list,root=0)
        sub_model_imagelist= numpy.array_split(model_imagelist, size)
        sub_model_imagelist=comm.scatter(sub_model_imagelist,root=0)
        if gcfcf is not None:
            sub_gcfcf = numpy.array_split(gcfcf,size)
            sub_gcfcf=comm.scatter(sub_gcfcf,root=0)
        isinstance(sub_vis_list[0], Visibility)
        image_results_list= [predict_list_serial_workflow(vis_list=[sub_vis_list[i]],
                     model_imagelist=[sub_model_imagelist[i]], vis_slices=vis_slices,
                     facets=facets, context=context, gcfcf=gcfcf, **kwargs)[0]
                for i, _ in enumerate(sub_vis_list)]
        #print(image_results_list)

        image_results_list=comm.gather(image_results_list,root=0)
        if rank == 0:
            #image_results_list_list=[x for x in image_results_list_list if x]
            image_results_list=numpy.concatenate(image_results_list)
        else:
            image_results_list=list()

    return image_results_list

def invert_list_mpi_workflow(vis_list, template_model_imagelist, context, dopsf=False, normalize=True,
                                    facets=1, vis_slices=1, gcfcf=None,
                                comm=MPI.COMM_WORLD, **kwargs):
    """ Sum results from invert, iterating over the scattered image and vis_list

    :param vis_list: Only full for rank==0
    :param template_model_imagelist: Model used to determine image parameters
    (in rank=0)
    :param dopsf: Make the PSF instead of the dirty image
    :param facets: Number of facets
    :param normalize: Normalize by sumwt
    :param vis_slices: Number of slices
    :param context: Imaging context
    :param gcfcg: tuple containing grid correction and convolution function (in
    rank=0)
    :param comm:MPI Communicator
    :param kwargs: Parameters for functions in components
    :return: List of (image, sumwt) tuple
   """

   # NOTE: Be careful with normalization as normalizing parts is not the
   # same as normalizing the whole, normalization happens for each image in a
   # frequency window (in this versio we only parallelize at freqwindows
    def concat_tuples(list_of_tuples):
        if len(list_of_tuples)<2:
            result_list=list_of_tuples
        else:
            result_list=list_of_tuples[0]
            for l in list_of_tuples[1:]:
                result_list+=l
        return result_list

    rank = comm.Get_rank()
    size = comm.Get_size()
    log.info('%d: In invert_list_mpi_workflow: %d elements in vis_list %d in model' %
          (rank,len(vis_list),len(template_model_imagelist)))

    assert get_parameter(kwargs, "use_serial_invert", True),"Only freq paralellization implemented"
    #if get_parameter(kwargs, "use_serial_invert", False):
    if get_parameter(kwargs, "use_serial_invert", True):
        from workflows.serial.imaging.imaging_serial import invert_list_serial_workflow

        results_vislist = list()
        # Distribute visibilities and model by freq
        sub_vis_list= numpy.array_split(vis_list, size)
        sub_vis_list=comm.scatter(sub_vis_list,root=0)
        sub_template_model_imagelist= numpy.array_split(template_model_imagelist, size)
        sub_template_model_imagelist=comm.scatter(sub_template_model_imagelist,root=0)
        if gcfcf is not None:
            sub_gcfcf = numpy.array_split(gcfcf,size)
            sub_gcfcf=comm.scatter(sub_gcfcf,root=0)
        isinstance(sub_vis_list[0], Visibility)
        sub_results_vislist = [invert_list_serial_workflow(vis_list=[sub_vis_list[i]],
                     template_model_imagelist=[sub_template_model_imagelist[i]],
                     context=context, dopsf=dopsf, normalize=normalize,
                                                       vis_slices=vis_slices,
                     facets=facets, gcfcf=gcfcf, **kwargs)[0]
                for i, _ in enumerate(sub_vis_list)]
        #print("%d sub_results_vislist" %rank,sub_results_vislist)
        results_vislist=comm.gather(sub_results_vislist,root=0)
        #print("%d results_vislist before concatenate"%rank,results_vislist)
        if rank == 0:
            #image_results_list_list=[x for x in image_results_list_list if x]
            #results_vislist=numpy.concatenate(results_vislist)
            # TODO: concatenate dos not concatenate well a list of tuples
            # it returns a 2d array instead of a concatenated list of tuples

            results_vislist =concat_tuples(results_vislist)
        else:
            results_vislist=list()
    #print("%d results_vislist"%rank,results_vislist)
    return results_vislist

def residual_list_mpi_workflow(vis, model_imagelist, context='2d', gcfcf=None,comm=MPI.COMM_WORLD, **kwargs):
    """ Create a graph to calculate residual image using w stacking and faceting

    :param context:
    :param vis: rank0
    :param model_imagelist: Model used to determine image parameters rank0
    :param gcfcg: tuple containing grid correction and convolution function
    :param kwargs: Parameters for functions in components
    :return:
    """
    rank = comm.Get_rank()
    size = comm.Get_size()
    model_vis = zero_list_mpi_workflow(vis)
    log.info('%d: In residual_list_mpi_workflow vis len %d model_imagelist len %d model_vis len %d'
          %(rank,len(vis),len(model_imagelist),len(model_vis)))

    model_vis = predict_list_mpi_workflow(model_vis, model_imagelist,
                                          context=context,gcfcf=gcfcf, **kwargs)
    residual_vis = subtract_list_mpi_workflow(vis, model_vis)
    return invert_list_mpi_workflow(residual_vis, model_imagelist, dopsf=False,
                                    normalize=True,
                                    context=context,gcfcf=gcfcf,
                                       **kwargs)

def restore_list_mpi_workflow(model_imagelist, psf_imagelist,
                              residual_imagelist,comm=MPI.COMM_WORLD, **kwargs):
    """ Create a graph to calculate the restored image

    :param model_imagelist: Model list (rank0)
    :param psf_imagelist: PSF list (rank0)
    :param residual_imagelist: Residual list (rank0)
    :param kwargs: Parameters for functions in components
    :return:
    """
    from workflows.serial.imaging.imaging_serial import restore_list_serial_workflow_nosumwt
    rank = comm.Get_rank()
    size = comm.Get_size()
    #TODO Parallelize! and check the dask version, it removes sumwt component
    # to reduce communication
    if residual_imagelist is None:
        residual_imagelist = []

    if rank==0:
        psf_list = remove_sumwt(psf_imagelist)
        if len(residual_imagelist) > 0:
            residual_list = remove_sumwt(residual_imagelist)
        else:
            residual_list = residual_imagelist
    else:
        psf_list=list()
        residual_list=list()

    sub_model_imagelist=numpy.array_split(model_imagelist,size)
    sub_model_imagelist=comm.scatter(sub_model_imagelist,root=0)
    sub_psf_list=numpy.array_split(psf_list,size)
    sub_psf_list=comm.scatter(sub_psf_list,root=0)
    sub_residual_list=numpy.array_split(residual_list,size)
    sub_residual_list=comm.scatter(sub_residual_list,root=0)

    sub_result_list=restore_list_serial_workflow_nosumwt(sub_model_imagelist,
                                         sub_psf_list,
                                         sub_residual_list)
    #sub_result_list=[restore_cube(sub_model_imagelist[i], sub_psf_list[i],
    #                     sub_residual_list[i], **kwargs)
    #            for i, _ in enumerate(sub_model_imagelist)]
    result_list=comm.gather(sub_result_list,root=0)
    if rank==0:
        # this is a list of tuples too, we may need to call my function
        result_list=numpy.concatenate(result_list)
    else:
        result_list=list()
    return result_list

def deconvolve_list_mpi_workflow(dirty_list, psf_list, model_imagelist,
                                        prefix='',
                                        mask=None,comm=MPI.COMM_WORLD, **kwargs):
    """Create a graph for deconvolution, adding to the model

    :param dirty_list: in rank0
    :param psf_list: in rank0
    :param model_imagelist: in rank0
    :param prefix: Informative prefix to log messages
    :param mask: Mask for deconvolution
    :param comm: MPI communicator
    :param kwargs: Parameters for functions in components
    :return: graph for the deconvolution
    """
    rank = comm.Get_rank()
    size = comm.Get_size()
    nchan = len(dirty_list)
    log.info('%d: deconvolve_list_mpi_workflow: dirty_list len %d psf_list len %d model_imagelist len %d' %(rank,len(dirty_list),len(psf_list),len(model_imagelist)))

    nmoment = get_parameter(kwargs, "nmoment", 0)
    assert get_parameter(kwargs, "use_serial_clean", True),"Only serial deconvolution implemented"
    if get_parameter(kwargs, "use_serial_clean", True):
        from workflows.serial.imaging.imaging_serial import deconvolve_list_serial_workflow
        if rank==0:
            assert isinstance(model_imagelist, list), model_imagelist
            result_list=deconvolve_list_serial_workflow (dirty_list, psf_list, model_imagelist, prefix=prefix, mask=mask, **kwargs)
        else:
            result_list=list()
    return result_list

def deconvolve_list_channel_mpi_workflow(dirty_list, psf_list, model_imagelist,
                                         subimages, comm=MPI.COMM_WORLD,**kwargs):
    """Create a graph for deconvolution by channels, adding to the model

    Does deconvolution channel by channel.
    :param subimages: MONTSE: number of subimages (= freqchannels?)
    :param dirty_list: in rank=0
    :param psf_list: Must be the size of a facet in rank=0
    :param model_imagelist: Current model in rank=0
    :param kwargs: Parameters for functions in components
    :return:
    """

    def deconvolve_subimage(dirty, psf):
        assert isinstance(dirty, Image)
        assert isinstance(psf, Image)
        comp = deconvolve_cube(dirty, psf, **kwargs)
        return comp[0]

    def add_model(sum_model, model):
        assert isinstance(output, Image)
        assert isinstance(model, Image)
        sum_model.data += model.data
        return sum_model

    rank = comm.Get_rank()
    size = comm.Get_size()
    if rank==0:
        output = create_empty_image_like(model_imagelist)
    dirty_lists = image_scatter_channels(dirty_list[0], subimages=subimages)
    sub_dirty_lists = numpy.array_split(dirty_lists,size)
    sub_dirty_lists = comm.scatter(sub_dirty_lists, root=0)
    psf_list_im=comm.Bcast(psf_list[0],root=0)

    sub_results = [deconvolve_subimage(dirty_list, psf_list_im)
               for dirty_list in sub_dirty_lists]
    results=comm.gather(sub_results,root=0)
    # NOTE: This is same as in invert, not scalable, we should use a reduction
    # instead but I don't understand image_gather_channels ...
    if rank==0:
        results=numpy.concatenate(results)
        result = image_gather_channels(results, output, subimages=subimages)
        result = add_model(result, model_imagelist)
    else:
        result=None
    return result

def weight_list_mpi_workflow(vis_list, model_imagelist, gcfcf=None,
                             weighting='uniform',comm=MPI.COMM_WORLD, **kwargs):
    """ Weight the visibility data

    This is done collectively so the weights are summed over all vis_lists and then
    corrected

    :param vis_list:
    :param model_imagelist: Model required to determine weighting parameters
    :param weighting: Type of weighting
    :param kwargs: Parameters for functions in graphs
    :return: List of vis_graphs
   """
    centre = len(model_imagelist) // 2
    rank = comm.Get_rank()
    size = comm.Get_size()

    if gcfcf is None:
        if rank==0:
            gcfcf = [create_pswf_convolutionfunction(model_imagelist[centre])]
    gcfcf = comm.bcast(gcfcf,root=0)

    def grid_wt(vis, model, g):
        if vis is not None:
            if model is not None:
                griddata = create_griddata_from_image(model)
                griddata = grid_weight_to_griddata(vis, griddata, g[0][1])
                return griddata
            else:
                return None
        else:
            return None

    sub_vis_list= numpy.array_split(vis_list, size)
    sub_vis_list=comm.scatter(sub_vis_list,root=0)
    sub_model_imagelist= numpy.array_split(model_imagelist, size)
    sub_model_imagelist=comm.scatter(sub_model_imagelist,root=0)

    sub_weight_list = [grid_wt(sub_vis_list[i], sub_model_imagelist[i], gcfcf)
                   for i in range(len(sub_vis_list))]

    weight_list=comm.allgather(sub_weight_list,root=0)
    weight_list=numpy.concatenate(weight_list)
    # This reduces a grid for each freqwin to a single grid  AllReduce but it
    # is not a simple reduction but I think it can be done with  4 reductions
    # Actually it seems I only use the  gd[0].data which would be the first
    # reduction only! TODO: try it out
    ##merged_sub_weight_grid = griddata_merge_weights(sub_weight_list)
    ##merged_weight_grid = comm.Allreduce(merged_sub_weight_grid[0].data,MPI.SUM)
    # Note that griddata_merge_weights does some averaging, if we need more
    # than the data component that would not be correct!
    merged_weight_grid = griddata_merge_weights(weight_list)

    def re_weight(vis, model, gd, g):
        if gd is not None:
            if vis is not None:
                # Ensure that the griddata has the right axes so that the convolution
                # function mapping works
                agd = create_griddata_from_image(model)
                agd.data = gd[0].data
                vis = griddata_reweight(vis, agd, g[0][1])
                return vis
            else:
                return None
        else:
            return vis

    sub_result = [re_weight(v, sub_model_imagelist[i], merged_weight_grid, gcfcf)
            for i, v in enumerate(sub_vis_list)]
    result=comm.gather(sub_result,root=0)
    if rank==0:
        result=numpy.concatenate(result)
    else:
        result=list()
    return result

def taper_list_mpi_workflow(vis_list, size_required,comm=MPI.COMM_WORLD ):
    """Taper to desired size

    :param vis_list:
    :param size_required:
    :return:
    """
    rank = comm.Get_rank()
    size = comm.Get_size()
    sub_vis_list= numpy.array_split(vis_list, size)
    sub_vis_list=comm.scatter(sub_vis_list,root=0)
    sub_result = [taper_visibility_gaussian(v, beam=size_required) for v in sub_vis_list]
    result=comm.gather(sub_result,root=0)
    if rank==0:
        result=numpy.concatenate(result)
    else:
        result=list()

    return result

def zero_list_mpi_workflow(vis_list,comm=MPI.COMM_WORLD ):
    """ Initialise vis to zero: creates new data holders

    :param vis_list:
    :return: List of vis_lists
   """

    def zero(vis):
        if vis is not None:
            zerovis = copy_visibility(vis)
            zerovis.data['vis'][...] = 0.0
            return zerovis
        else:
            return None

    rank = comm.Get_rank()
    size = comm.Get_size()
    sub_vis_list= numpy.array_split(vis_list, size)
    sub_vis_list=comm.scatter(sub_vis_list,root=0)
    sub_result = [zero(v) for v in sub_vis_list]
    result=comm.gather(sub_result,root=0)
    if rank==0:
        result=numpy.concatenate(result)
    else:
        result=list()
    return result

def subtract_list_mpi_workflow(vis_list, model_vislist,comm=MPI.COMM_WORLD ):
    """ Initialise vis to zero

    :param vis_list:
    :param model_vislist: Model to be subtracted
    :return: List of vis_lists
   """

    def subtract_vis(vis, model_vis):
        if vis is not None and model_vis is not None:
            assert vis.vis.shape == model_vis.vis.shape,vis
            subvis = copy_visibility(vis)
            subvis.data['vis'][...] -= model_vis.data['vis'][...]
            return subvis
        else:
            return None
    rank = comm.Get_rank()
    size = comm.Get_size()
    print("%d: In substract: vis_list" %rank,vis_list)
    sub_vis_list= numpy.array_split(vis_list, size)
    sub_vis_list=comm.scatter(sub_vis_list,root=0)
    sub_model_vislist= numpy.array_split(model_vislist, size)
    sub_model_vislist=comm.scatter(sub_model_vislist,root=0)
    print("%d In substract" %rank,"model:",sub_model_vislist,"vis:",sub_vis_list)
    sub_result = [subtract_vis(vis=sub_vis_list[i], model_vis=sub_model_vislist[i])
            for i in range(len(sub_vis_list))]
    result=comm.gather(sub_result,root=0)
    if rank==0:
        result=numpy.concatenate(result)
    else:
        result=list()
    return result


text:
import numpy as np
import torch
import torch.utils.data
from torch.autograd import Variable
import torch.nn.functional as F
import ae as models
from mmd import mix_rbf_mmd2
import math
import time

# sigma for MMD
base = 1.0
sigma_list = [1, 2, 4, 8, 16]
sigma_list = [sigma / base for sigma in sigma_list]

from imblearn.over_sampling import RandomOverSampler
imblearn_seed = 0

def training(dataset_list, cluster_pairs, nn_paras):
    """ Training an autoencoder to remove batch effects
    Args:
        dataset_list: list of datasets for batch correction
        cluster_pairs: pairs of similar clusters with weights
        nn_paras: parameters for neural network training
    Returns:
        model: trained autoencoder
        loss_total_list: list of total loss
        loss_reconstruct_list: list of reconstruction loss
        loss_transfer_list: list of transfer loss
    """
    # load nn parameters
    batch_size = nn_paras['batch_size']
    num_epochs = nn_paras['num_epochs']
    num_inputs = nn_paras['num_inputs']
    code_dim = nn_paras['code_dim']
    cuda = nn_paras['cuda']

    # training data for autoencoder, construct a DataLoader for each cluster
    cluster_loader_dict = {}
    for i in range(len(dataset_list)):
        gene_exp = dataset_list[i]['gene_exp'].transpose()
        cluster_labels = dataset_list[i]['cluster_labels']  # cluster labels do not overlap between datasets
        unique_labels = np.unique(cluster_labels)
        # Random oversampling based on cell cluster sizes
        gene_exp, cluster_labels = RandomOverSampler(random_state=imblearn_seed).fit_sample(gene_exp, cluster_labels)

        # construct DataLoader list
        for j in range(len(unique_labels)):
            idx = cluster_labels == unique_labels[j]
            if cuda:
                torch_dataset = torch.utils.data.TensorDataset(
                    torch.FloatTensor(gene_exp[idx,:]).cuda(), torch.LongTensor(cluster_labels[idx]).cuda())
            else:
                torch_dataset = torch.utils.data.TensorDataset(
                    torch.FloatTensor(gene_exp[idx, :]), torch.LongTensor(cluster_labels[idx]))
            data_loader = torch.utils.data.DataLoader(torch_dataset, batch_size=batch_size,
                                                        shuffle=True, drop_last=True)
            cluster_loader_dict[unique_labels[j]] = data_loader

    # create model
    if code_dim == 20:
        model = models.autoencoder_20(num_inputs=num_inputs)
    elif code_dim == 2:
        model = models.autoencoder_2(num_inputs=num_inputs)
    else:
        model = models.autoencoder_20(num_inputs=num_inputs)
    if cuda:
        model.cuda()

    # training
    loss_total_list = []  # list of total loss
    loss_reconstruct_list = []
    loss_transfer_list = []
    for epoch in range(1, num_epochs + 1):
        avg_loss, avg_reco_loss, avg_tran_loss = training_epoch(epoch, model, cluster_loader_dict, cluster_pairs, nn_paras)
        # terminate early if loss is nan
        if math.isnan(avg_reco_loss) or math.isnan(avg_tran_loss):
            return [], model, [], [], []
        loss_total_list.append(avg_loss)
        loss_reconstruct_list.append(avg_reco_loss)
        loss_transfer_list.append(avg_tran_loss)

    return model, loss_total_list, loss_reconstruct_list, loss_transfer_list

def training_epoch(epoch, model, cluster_loader_dict, cluster_pairs, nn_paras):
    """ Training an epoch
        Args:
            epoch: number of the current epoch
            model: autoencoder
            cluster_loader_dict: dict of DataLoaders indexed by clusters
            cluster_pairs: pairs of similar clusters with weights
            nn_paras: parameters for neural network training
        Returns:
            avg_total_loss: average total loss of mini-batches
            avg_reco_loss: average reconstruction loss of mini-batches
            avg_tran_loss: average transfer loss of mini-batches
        """
    log_interval = nn_paras['log_interval']
    # load nn parameters
    base_lr = nn_paras['base_lr']
    lr_step = nn_paras['lr_step']
    num_epochs = nn_paras['num_epochs']
    l2_decay = nn_paras['l2_decay']
    gamma = nn_paras['gamma']
    cuda = nn_paras['cuda']

    # step decay of learning rate
    learning_rate = base_lr / math.pow(2, math.floor(epoch / lr_step))
    # regularization parameterbetween two losses
    gamma_rate = 2 / (1 + math.exp(-10 * (epoch) / num_epochs)) - 1
    gamma = gamma_rate * gamma

    if epoch % log_interval == 0:
        print('{:}, Epoch {}, learning rate {:.3E}, gamma {:.3E}'.format(
                time.asctime(time.localtime()), epoch, learning_rate, gamma))

    optimizer = torch.optim.Adam([
        {'params': model.encoder.parameters()},
        {'params': model.decoder.parameters()},
    ], lr=learning_rate, weight_decay=l2_decay)

    model.train()

    iter_data_dict = {}
    for cls in cluster_loader_dict:
        iter_data = iter(cluster_loader_dict[cls])
        iter_data_dict[cls] = iter_data
    # use the largest dataset to define an epoch
    num_iter = 0
    for cls in cluster_loader_dict:
        num_iter = max(num_iter, len(cluster_loader_dict[cls]))

    total_loss = 0
    total_reco_loss = 0
    total_tran_loss = 0
    num_batches = 0

    for it in range(0, num_iter):
        data_dict = {}
        label_dict = {}
        code_dict = {}
        reconstruct_dict = {}
        for cls in iter_data_dict:
            data, labels = iter_data_dict[cls].next()
            data_dict[cls] = data
            label_dict[cls] = labels
            if it % len(cluster_loader_dict[cls]) == 0:
                iter_data_dict[cls] = iter(cluster_loader_dict[cls])
            data_dict[cls] = Variable(data_dict[cls])
            label_dict[cls] = Variable(label_dict[cls])

        for cls in data_dict:
            code, reconstruct = model(data_dict[cls])
            code_dict[cls] = code
            reconstruct_dict[cls] = reconstruct

        optimizer.zero_grad()

        # transfer loss for cluster pairs in cluster_pairs matrix
        loss_transfer = torch.FloatTensor([0])
        if cuda:
            loss_transfer = loss_transfer.cuda()
        for i in range(cluster_pairs.shape[0]):
            cls_1 = int(cluster_pairs[i,0])
            cls_2 = int(cluster_pairs[i,1])
            if cls_1 not in code_dict or cls_2 not in code_dict:
                continue
            mmd2_D = mix_rbf_mmd2(code_dict[cls_1], code_dict[cls_2], sigma_list)
            loss_transfer += mmd2_D * cluster_pairs[i,2]

        # reconstruction loss for all clusters
        loss_reconstruct = torch.FloatTensor([0])
        if cuda:
            loss_reconstruct = loss_reconstruct.cuda()
        for cls in data_dict:
            loss_reconstruct += F.mse_loss(reconstruct_dict[cls], data_dict[cls])

        loss = loss_reconstruct + gamma * loss_transfer

        loss.backward()
        optimizer.step()

        # update total loss
        num_batches += 1
        total_loss += loss.data.item()
        total_reco_loss += loss_reconstruct.data.item()
        total_tran_loss += loss_transfer.data.item()

    avg_total_loss = total_loss / num_batches
    avg_reco_loss = total_reco_loss / num_batches
    avg_tran_loss = total_tran_loss / num_batches

    if epoch % log_interval == 0:
        print('Avg_loss {:.3E}\t Avg_reconstruct_loss {:.3E}\t Avg_transfer_loss {:.3E}'.format(
            avg_total_loss, avg_reco_loss, avg_tran_loss))
    return avg_total_loss, avg_reco_loss, avg_tran_loss

def testing(model, dataset_list, nn_paras):
    """ Training an epoch
    Args:
        model: autoencoder
        dataset_list: list of datasets for batch correction
        nn_paras: parameters for neural network training
    Returns:
        code_list: list pf embedded codes
    """

    # load nn parameters
    cuda = nn_paras['cuda']

    data_loader_list = []
    num_cells = []
    for dataset in dataset_list:
        torch_dataset = torch.utils.data.TensorDataset(
            torch.FloatTensor(dataset['gene_exp'].transpose()), torch.LongTensor(dataset['cell_labels']))
        data_loader = torch.utils.data.DataLoader(torch_dataset, batch_size=len(dataset['cell_labels']),
                                                    shuffle=False)
        data_loader_list.append(data_loader)
        num_cells.append(len(dataset["cell_labels"]))

    model.eval()

    code_list = [] # list pf embedded codes
    for i in range(len(data_loader_list)):
        idx = 0
        with torch.no_grad():
            for data, labels in data_loader_list[i]:
                if cuda:
                    data, labels = data.cuda(), labels.cuda()
                code_tmp, _ = model(data)
                code_tmp = code_tmp.cpu().numpy()
                if idx == 0:
                    code = np.zeros((code_tmp.shape[1], num_cells[i]))
                code[:, idx:idx + code_tmp.shape[0]] = code_tmp.T
                idx += code_tmp.shape[0]
        code_list.append(code)

    return code_list

ext: py

reason: import numpy


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

# Reference: https://www.kaggle.com/rishabhiitbhu/eda-understanding-the-dataset-with-3d-plots

import os
import pdb
import cv2
import pandas as pd
import numpy as np
from PIL import Image
from pathlib import Path
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import axes3d, Axes3D

from lyft_dataset_sdk.lyftdataset import LyftDataset, LyftDatasetExplorer, Quaternion, view_points
from lyft_dataset_sdk.utils.data_classes import LidarPointCloud

from moviepy.editor import ImageSequenceClip
from tqdm import tqdm_notebook as tqdm

def render_scene(index):
    my_scene = lyftdata.scene[index]
    my_sample_token = my_scene["first_sample_token"]
    lyftdata.render_sample(my_sample_token)

def display_frame_statistics(lidar_points, points=0.2):
    """
    Displays statistics for a single frame. Draws 3D plot of the lidar point cloud data and point cloud
    projections to various planes.

    Parameters
    ----------
    lidar_points: lidar data points
    points          : Fraction of lidar points to use. Defaults to `0.2`, e.g. 20%.
    """

    points_step = int(1. / points)
    point_size = 0.01 * (1. / points)
    pc_range = range(0, lidar_points.shape[1], points_step)
    pc_frame = lidar_points[:, pc_range]
    def draw_point_cloud(ax, title, axes=[0, 1, 2]):
        """Convenient method for drawing various point cloud projections as a part of frame statistics"""
        ax.set_facecolor('black')
        ax.grid(False)
        ax.scatter(*pc_frame[axes, :], s=point_size, c='white', cmap='grey')
        if len(axes) == 3: # 3D configs
            text_color = 'white'
            ax.set_xlim3d([-10, 30])
            ax.set_ylim3d(*axes_limits[axes[1]])
            ax.set_zlim3d(*axes_limits[axes[2]])
            ax.set_zlabel('{} axis'.format(axes_str[axes[2]]), color='white')
        else: # 2D configs
            text_color = 'black' # the `figure` is white
            ax.set_xlim(*axes_limits[axes[0]])
            ax.set_ylim(*axes_limits[axes[1]])
        ax.set_title(title, color=text_color)
        ax.set_xlabel('{} axis'.format(axes_str[axes[0]]), color=text_color)
        ax.set_ylabel('{} axis'.format(axes_str[axes[1]]), color=text_color)

    # Draw point cloud data as 3D plot
    f2 = plt.figure(figsize=(15, 8))
    ax2 = f2.add_subplot(111, projection='3d')
    # make the panes transparent
    ax2.xaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    ax2.yaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    ax2.zaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    draw_point_cloud(ax2, '3D plot')

    plt.show()
    # Draw point cloud data as plane projections
    f, ax3 = plt.subplots(3, 1, figsize=(15, 25))
#     f.set_facecolor('black')
    draw_point_cloud(
        ax3[0],
        'XZ projection (Y = 0)', #, the car is moving in direction left to right', ?
        axes=[0, 2] # X and Z axes
    )
    draw_point_cloud(
        ax3[1],
        'XY projection (Z = 0)', #, the car is moving in direction left to right',?
        axes=[0, 1] # X and Y axes
    )
    draw_point_cloud(
        ax3[2],
        'YZ projection (X = 0)', #, the car is moving towards the graph plane', ?
        axes=[1, 2] # Y and Z axes
    )
    plt.show()

def get_lidar_points(lidar_token):
    '''Get lidar point cloud in the frame of the ego vehicle'''
    sd_record = lyftdata.get("sample_data", lidar_token)
    sensor_modality = sd_record["sensor_modality"]

    # Get aggregated point cloud in lidar frame.
    sample_rec = lyftdata.get("sample", sd_record["sample_token"])
    chan = sd_record["channel"]
    ref_chan = "LIDAR_TOP"
    pc, times = LidarPointCloud.from_file_multisweep(
        lyftdata, sample_rec, chan, ref_chan, num_sweeps=1
    )
    # Compute transformation matrices for lidar point cloud
    cs_record = lyftdata.get("calibrated_sensor", sd_record["calibrated_sensor_token"])
    pose_record = lyftdata.get("ego_pose", sd_record["ego_pose_token"])
    vehicle_from_sensor = np.eye(4)
    vehicle_from_sensor[:3, :3] = Quaternion(cs_record["rotation"]).rotation_matrix
    vehicle_from_sensor[:3, 3] = cs_record["translation"]

    ego_yaw = Quaternion(pose_record["rotation"]).yaw_pitch_roll[0]
    rot_vehicle_flat_from_vehicle = np.dot(
        Quaternion(scalar=np.cos(ego_yaw / 2), vector=[0, 0, np.sin(ego_yaw / 2)]).rotation_matrix,
        Quaternion(pose_record["rotation"]).inverse.rotation_matrix,
    )
    vehicle_flat_from_vehicle = np.eye(4)
    vehicle_flat_from_vehicle[:3, :3] = rot_vehicle_flat_from_vehicle
    points = view_points(
        pc.points[:3, :], np.dot(vehicle_flat_from_vehicle, vehicle_from_sensor), normalize=False
    )
    return points

def plot_box(box, axis, view, colors, normalize=False, linewidth=1.0):
    '''Plot boxes in the 3d figure'''
    corners = view_points(box.corners(), view, normalize=normalize)#
    def draw_rect(selected_corners, color):
        prev = selected_corners[-1]
        for corner in selected_corners:
            axis.plot([prev[0], corner[0]], [prev[1], corner[1]], [prev[2], corner[2]], color=color, linewidth=linewidth)
            prev = corner

    # Draw the sides
    for i in range(4):
        axis.plot(
            [corners.T[i][0], corners.T[i + 4][0]],
            [corners.T[i][1], corners.T[i + 4][1]],
            [corners.T[i][2], corners.T[i + 4][2]],
            color=colors[2],
            linewidth=linewidth,
        )
    # Draw front (first 4 corners) and rear (last 4 corners) rectangles(3d)/lines(2d)
    draw_rect(corners.T[:4], colors[0]) #4x3
    draw_rect(corners.T[4:], colors[1])

def draw_3d_plot(idx, lidar_token):
    '''Plot the lidar + annotations on a 3D figure'''
    # sample lidar point cloud
    lidar_points = get_lidar_points(lidar_token)
    points = 0.5 # fraction of lidar_points to plot, to reduce the clutter
    points_step = int(1. / points)
    pc_range = range(0, lidar_points.shape[1], points_step)
    lidar_points = lidar_points[:, pc_range]

    # Get boxes, instead of current sensor's coordinate frame, use vehicle frame which is aligned to z-plane in world
    _, boxes, _ = lyftdata.get_sample_data(
        lidar_token, flat_vehicle_coordinates=True
    )
    fig = plt.figure(figsize=(15, 8))
    ax = fig.add_subplot(111, projection='3d')
    point_size = 0.01 * (1. / points) # size of the dots on plot
    ax.set_facecolor('black')
    ax.grid(False)
    ax.scatter(*lidar_points, s=point_size, c='white', cmap='gray')
    for box in boxes:
        c = np.array(lyftdata.explorer.get_color(box.name)) / 255.0
        plot_box(box, ax, view=np.eye(3), colors=(c, c, c), linewidth=0.5)
    ax.set_xlim3d(-40, 40)
    ax.set_ylim3d(-40, 40)
    ax.set_zlim3d(-4, 40)

    # make the panes transparent
    ax.xaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    ax.yaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    ax.zaxis.set_pane_color((1.0, 1.0, 1.0, 0.0))
    ax.set_title(lidar_token, color='white')
    filename = DATA_PATH+r'tmp\frame_{0:0>4}.png'.format(idx)
    plt.savefig(filename)
    plt.close()
    return filename

## load data, you only have to do it once

DATA_PATH =r'C:\Users\storm\Documents\3d-object-detection-for-autonomous-vehicles\train_data'+os.sep
train = pd.read_csv(DATA_PATH + 'train.csv')
sample_submission = pd.read_csv(DATA_PATH + 'sample_submission.csv')
lyftdata = LyftDataset(data_path=DATA_PATH, json_path=DATA_PATH+os.sep+'data')

## show point cloud interactively

cat_token = lyftdata.category[0]['token']
train = pd.read_csv(DATA_PATH+os.sep+'train.csv')
token0 = train.iloc[0]['Id']
my_sample = lyftdata.get('sample', token0)
lyftdata.render_sample_3d_interactive(my_sample['token'], render_sample=False)

## sensor information

sensor = 'CAM_FRONT'
cam_front = lyftdata.get('sample_data', my_sample['data'][sensor])
img = Image.open(DATA_PATH+cam_front['filename'])
lyftdata.render_sample_data(cam_front['token'], with_anns=False)

## look at the LIDAR data associated with my_sample

lidar_top = lyftdata.get('sample_data', my_sample['data']['LIDAR_TOP']) # selecting LIDAR_TOP out of all LIDARs
pc = LidarPointCloud.from_file(Path(DATA_PATH+lidar_top['filename']))

## then visualize lidar points

axes_limits = [
    [-30, 50], # X axis range
    [-30, 20], # Y axis range
    [-3, 10]   # Z axis range
]
axes_str = ['X', 'Y', 'Z']
display_frame_statistics(pc.points, points=0.5)
lyftdata.render_sample_data(lidar_top['token'], with_anns=False)

## Ego pose

lyftdata.get('ego_pose', cam_front['ego_pose_token'])

## caliberated_sensor

lyftdata.get('calibrated_sensor', cam_front['calibrated_sensor_token'])

## Annotations

my_annotation = lyftdata.get('sample_annotation', my_sample['anns'][0])
my_box = lyftdata.get_box(my_annotation['token'])
lyftdata.render_annotation(my_annotation['token'], margin=10)

## Attributes

my_attribute1 = lyftdata.get('attribute', my_annotation['attribute_tokens'][0])
my_attribute2 = lyftdata.get('attribute', my_annotation['attribute_tokens'][1])

## Instances

my_instance = lyftdata.get('instance', my_annotation['instance_token'])
lyftdata.render_instance(my_instance['token'])
print("First annotated sample of this instance:")
lyftdata.render_annotation(my_instance['first_annotation_token'])
print("Last annotated sample of this instance")
lyftdata.render_annotation(my_instance['last_annotation_token'])
lyftdata.render_sample(token0)

## 3D visualization of a scene
## Be sure to make a folder named tmp in data_path

my_scene = lyftdata.get('scene',  my_sample['scene_token'])
first_sample_token = my_scene['first_sample_token']
sample = lyftdata.get('sample', first_sample_token)
lidar_token = sample['data']['LIDAR_TOP']
filename = draw_3d_plot(0, lidar_token)
Image.open(filename)


text:
r"""
Hopf algebras with basis
"""
# ****************************************************************************
#  Copyright (C) 2008 Teresa Gomez-Diaz (CNRS) <Teresa.Gomez-Diaz@univ-mlv.fr>
#  Copyright (C) 2008-2011 Nicolas M. Thiery <nthiery at users.sf.net>
#
#  Distributed under the terms of the GNU General Public License (GPL)
#                  https://www.gnu.org/licenses/
# *****************************************************************************

from sage.categories.category_with_axiom import CategoryWithAxiom_over_base_ring
from sage.categories.tensor import TensorProductsCategory
from sage.misc.abstract_method import abstract_method
from sage.misc.cachefunc import cached_method
from sage.misc.lazy_attribute import lazy_attribute
from sage.misc.lazy_import import LazyImport

class HopfAlgebrasWithBasis(CategoryWithAxiom_over_base_ring):
    """
    The category of Hopf algebras with a distinguished basis

    EXAMPLES::

        sage: C = HopfAlgebrasWithBasis(QQ)
        sage: C
        Category of hopf algebras with basis over Rational Field
        sage: C.super_categories()
        [Category of hopf algebras over Rational Field,
         Category of bialgebras with basis over Rational Field]

    We now show how to use a simple Hopf algebra, namely the group algebra of the dihedral group
    (see also AlgebrasWithBasis)::

        sage: A = C.example(); A
        An example of Hopf algebra with basis: the group algebra of the Dihedral group of order 6 as a permutation group over Rational Field
        sage: A.__custom_name = "A"
        sage: A.category()
        Category of finite dimensional hopf algebras with basis over Rational Field

        sage: A.one_basis()
        ()
        sage: A.one()
        B[()]

        sage: A.base_ring()
        Rational Field
        sage: A.basis().keys()
        Dihedral group of order 6 as a permutation group

        sage: [a,b] = A.algebra_generators()
        sage: a, b
        (B[(1,2,3)], B[(1,3)])
        sage: a^3, b^2
        (B[()], B[()])
        sage: a*b
        B[(1,2)]

        sage: A.product           # todo: not quite ...
        <bound method MagmaticAlgebras.WithBasis.ParentMethods._product_from_product_on_basis_multiply of A>
        sage: A.product(b,b)
        B[()]

        sage: A.zero().coproduct()
        0
        sage: A.zero().coproduct().parent()
        A # A
        sage: a.coproduct()
        B[(1,2,3)] # B[(1,2,3)]

        sage: TestSuite(A).run(verbose=True)
        running ._test_additive_associativity() . . . pass
        running ._test_an_element() . . . pass
        running ._test_antipode() . . . pass
        running ._test_associativity() . . . pass
        running ._test_cardinality() . . . pass
        running ._test_category() . . . pass
        running ._test_characteristic() . . . pass
        running ._test_construction() . . . pass
        running ._test_distributivity() . . . pass
        running ._test_elements() . . .
          Running the test suite of self.an_element()
          running ._test_category() . . . pass
          running ._test_eq() . . . pass
          running ._test_new() . . . pass
          running ._test_nonzero_equal() . . . pass
          running ._test_not_implemented_methods() . . . pass
          running ._test_pickling() . . . pass
          pass
        running ._test_elements_eq_reflexive() . . . pass
        running ._test_elements_eq_symmetric() . . . pass
        running ._test_elements_eq_transitive() . . . pass
        running ._test_elements_neq() . . . pass
        running ._test_eq() . . . pass
        running ._test_new() . . . pass
        running ._test_not_implemented_methods() . . . pass
        running ._test_one() . . . pass
        running ._test_pickling() . . . pass
        running ._test_prod() . . . pass
        running ._test_some_elements() . . . pass
        running ._test_zero() . . . pass
        sage: A.__class__
        <class 'sage.categories.examples.hopf_algebras_with_basis.MyGroupAlgebra_with_category'>
        sage: A.element_class
        <class 'sage.categories.examples.hopf_algebras_with_basis.MyGroupAlgebra_with_category.element_class'>

    Let us look at the code for implementing A::

        sage: A??                       # todo: not implemented

    TESTS::

        sage: TestSuite(A).run()
        sage: TestSuite(C).run()
    """

    def example(self, G = None):
        """
        Returns an example of algebra with basis::

            sage: HopfAlgebrasWithBasis(QQ['x']).example()
            An example of Hopf algebra with basis: the group algebra of the Dihedral group of order 6 as a permutation group over Univariate Polynomial Ring in x over Rational Field

        An other group can be specified as optional argument::

            sage: HopfAlgebrasWithBasis(QQ).example(SymmetricGroup(4))
            An example of Hopf algebra with basis: the group algebra of the Symmetric group of order 4! as a permutation group over Rational Field
        """
        from sage.categories.examples.hopf_algebras_with_basis import MyGroupAlgebra
        from sage.groups.perm_gps.permgroup_named import DihedralGroup
        if G is None:
            G = DihedralGroup(3)
        return MyGroupAlgebra(self.base_ring(), G)

#     This is only correct in the finite dimensional / graded case
#     def dual(self):
#         """
#         Returns the dual category

#         EXAMPLES:

#         The category of Hopf algebras over any field is self dual::

#             sage: C = HopfAlgebrasWithBasis(QQ)
#             sage: C.dual()
#             Category of hopf algebras with basis over Rational Field
#         """
#         return self

    FiniteDimensional = LazyImport('sage.categories.finite_dimensional_hopf_algebras_with_basis',
                                   'FiniteDimensionalHopfAlgebrasWithBasis')
    Filtered = LazyImport('sage.categories.filtered_hopf_algebras_with_basis',
                          'FilteredHopfAlgebrasWithBasis')
    Graded = LazyImport('sage.categories.graded_hopf_algebras_with_basis',
                        'GradedHopfAlgebrasWithBasis')
    Super = LazyImport('sage.categories.super_hopf_algebras_with_basis',
                       'SuperHopfAlgebrasWithBasis')

    class ParentMethods:

        @abstract_method(optional=True)
        def antipode_on_basis(self, x):
            """
            The antipode of the Hopf algebra on the basis (optional)

            INPUT:

             - ``x`` -- an index of an element of the basis of ``self``

            Returns the antipode of the basis element indexed by ``x``.

            If this method is implemented, then :meth:`antipode` is defined
            from this by linearity.

            EXAMPLES::

                sage: A = HopfAlgebrasWithBasis(QQ).example()
                sage: W = A.basis().keys(); W
                Dihedral group of order 6 as a permutation group
                sage: w = W.gen(0); w
                (1,2,3)
                sage: A.antipode_on_basis(w)
                B[(1,3,2)]
            """

        @lazy_attribute
        def antipode(self):
            """
            The antipode of this Hopf algebra.

            If :meth:`.antipode_basis` is available, this constructs the
            antipode morphism from ``self`` to ``self`` by extending it by
            linearity. Otherwise, :meth:`self.antipode_by_coercion` is used, if
            available.

            EXAMPLES::

                sage: A = HopfAlgebrasWithBasis(ZZ).example(); A
                An example of Hopf algebra with basis: the group algebra of the Dihedral group of order 6 as a permutation group over Integer Ring
                sage: A = HopfAlgebrasWithBasis(QQ).example()
                sage: [a,b] = A.algebra_generators()
                sage: a, A.antipode(a)
                (B[(1,2,3)], B[(1,3,2)])
                sage: b, A.antipode(b)
                (B[(1,3)], B[(1,3)])

            TESTS::

                sage: all(A.antipode(x) * x == A.one() for x in A.basis())
                True
            """
            if self.antipode_on_basis is not NotImplemented:
                # Should give the information that this is an anti-morphism of algebra
                return self._module_morphism(self.antipode_on_basis, codomain = self)
            elif hasattr(self, "antipode_by_coercion"):
                return self.antipode_by_coercion

        def _test_antipode(self, **options):
            r"""
            Test the antipode.

            An *antipode* `S` of a Hopf algebra is a linear endomorphism of the
            Hopf algebra that satisfies the following conditions (see
            :wikipedia:`HopfAlgebra`).

            - If `\mu` and `\Delta` denote the product and coproduct of the
              Hopf algebra, respectively, then `S` satisfies

              .. MATH::

                  \mu \circ (S \tensor 1) \circ \Delta = unit \circ counit
                  \mu \circ (1 \tensor S) \circ \Delta = unit \circ counit

            - `S` is an *anti*-homomorphism

            These properties are tested on :meth:`some_elements`.

            TESTS::

                sage: R = NonCommutativeSymmetricFunctions(QQ).ribbon()
                sage: R._test_antipode()

            ::

                sage: s = SymmetricFunctions(QQ).schur()
                sage: s._test_antipode()

            """
            tester = self._tester(**options)

            S = self.antipode

            IS = lambda x: self.sum(c * self.monomial(t1) * S(self.monomial(t2))
                                for ((t1, t2), c) in x.coproduct())

            SI = lambda x: self.sum(c * S(self.monomial(t1)) * self.monomial(t2)
                                for ((t1, t2), c) in x.coproduct())

            for x in tester.some_elements():

                # antipode is an anti-homomorphism
                for y in tester.some_elements():
                    tester.assertEqual(S(x) * S(y), S(y * x))

                # mu * (S # I) * delta == counit * unit
                tester.assertEqual(SI(x), self.counit(x) * self.one())

                # mu * (I # S) * delta == counit * unit
                tester.assertEqual(IS(x), self.counit(x) * self.one())

    class ElementMethods:
        pass

    class TensorProducts(TensorProductsCategory):
        """
        The category of hopf algebras with basis constructed by tensor product of hopf algebras with basis
        """

        @cached_method
        def extra_super_categories(self):
            """
            EXAMPLES::

                sage: C = HopfAlgebrasWithBasis(QQ).TensorProducts()
                sage: C.extra_super_categories()
                [Category of hopf algebras with basis over Rational Field]
                sage: sorted(C.super_categories(), key=str)
                [Category of hopf algebras with basis over Rational Field,
                 Category of tensor products of algebras with basis over Rational Field,
                 Category of tensor products of hopf algebras over Rational Field]
            """
            return [self.base_category()]

        class ParentMethods:
            # todo: antipode
            pass

        class ElementMethods:
            pass

max_stars_repo_name: LaisRast/sage

max_issues_repo_name: LaisRast/sage

max_forks_repo_name: LaisRast/sage

reason: from sage


text:
from numpy import product
from rest_framework import serializers
from product.serializers import ProductSerializer

from order.models import Order
from product.models import Product

class OrderSerializer(serializers.ModelSerializer):
    products = ProductSerializer(read_only=True, many=True)
    products_ids = serializers.PrimaryKeyRelatedField(
        queryset=Product.objects.all(), write_only=True, many=True)
    total = serializers.SerializerMethodField()

    def get_total(self, instance):
        total = sum([product.price for product in instance.products.all()])
        return total

    class Meta:
        model = Order
        fields = ["products", "total", 'products_ids', 'user']
        extra_kwargs = {'product': {'required': False}}

    def create(self, validated_data):
        product_data = validated_data.pop('products_ids')
        order = Order.objects.create(**validated_data)
        for product in product_data:
            order.products.add(product)
        return order

max_stars_repo_name: DanLivassan/bookstore

max_issues_repo_name: DanLivassan/bookstore

max_forks_repo_name: DanLivassan/bookstore

reason: from numpy


text:
import numpy as np
import chainer
from scipy.spatial.transform import Rotation

import chainer_graphics.transform as T

TESTCASES = 10

def test_euler2rot():
    angles = np.random.randn(TESTCASES, 3)
    R1 = Rotation.from_euler('xyz', angles).as_dcm()
    R2 = T.euler2rot(chainer.Variable(angles)).data
    assert(np.allclose(R1, R2))

def test_rot2euler():
    angles1 = np.random.randn(TESTCASES, 3)
    R1 = Rotation.from_euler('xyz', angles1).as_dcm()
    angles2 = T.rot2euler(chainer.Variable(R1)).data
    R2 = Rotation.from_euler('xyz', angles2).as_dcm()
    assert(np.allclose(R1, R2))

max_stars_repo_name: Idein/chainer-graphics

max_issues_repo_name: Idein/chainer-graphics

max_forks_repo_name: Idein/chainer-graphics

reason: import numpy,from scipy


text:
import numpy as np
from matplotlib.path import Path
from numpy import vstack, hstack

# Matlab原始版 inpolygon函數
def inpolygon_Matlab(x_point, y_point, x_area, y_area):
    """
    Reimplement inpolygon in matlab
    :type x_point: np.ndarray
    :type y_point: np.ndarray
    :type x_area: np.ndarray
    :type y_area: np.ndarray
    """
    # 合併xv和yv為頂點數組
    vertices = vstack((x_area, y_area)).T
    # 定義Path對象
    path = Path(vertices)
    # 把xq和yq合併為test_points
    test_points = hstack([x_point.reshape(x_point.size, -1), y_point.reshape(y_point.size, -1)])
    # 得到一個test_points是否嚴格在path內的mask，是bool值數組
    _in = path.contains_points(test_points)
    # 得到一個test_points是否在path內部或者在路徑上的mask
    _in_on = path.contains_points(test_points, radius=-1e-10)
    # 得到一個test_points是否在path路徑上的mask
    _on = _in ^ _in_on
    return _in_on, _on

# 簡化inpolygon函數
def inpolygon(xq, yq, xv, yv):
    """
    reimplement inpolygon in matlab
    :type xq: np.ndarray
    :type yq: np.ndarray
    :type xv: np.ndarray
    :type yv: np.ndarray
    """
    # 合併xv和yv為頂點數組
    vertices = np.vstack((xv, yv)).T
    # 定義Path對象
    path = Path(vertices)
    # 把xq和yq合併為test_points
    test_points = np.hstack([xq.reshape(xq.size, -1), yq.reshape(yq.size, -1)])
    # 得到一個test_points是否嚴格在path內的mask，是bool值數組
    _in = path.contains_points(test_points)
    # 得到一個test_points是否在path內部或者在路徑上的mask
    _in_on = path.contains_points(test_points, radius=-1e-10)

    return _in_on

if __name__ == '__main__':
    # inpolygon 函數用法
    xv = np.array([-4, 0, 4, 0])
    yv = np.array([0, 4, 0, -4])
    X = np.array([0, 1, 3.5, 4, 5])
    Y = np.array([0, 1, 0, 0, 0])

    _in, _on = inpolygon_Matlab(X, Y, xv, yv)

    print(_in)
    print(_on)

reason: import numpy,from numpy


text:
# -*- coding: utf-8 -*-
"""
Created on Fri Mar 10 11:15:04 2017

@author: nchitta
"""
###Change Log 3/13/2017 Changes have been made to get optimal 4411B settings to be###
###able to get better match to the output using some statistical techniques.      ###
#####################################################################################
#   This code is trying to control all the instruments for basic GPIB programming   #
#   to conduct the random error experiment.                                         #
#####################################################################################
# Refer randerr experiment.txt in /home/chitta/Desktop/Python/InstCon/ for backgrnd #
#####################################################################################
# Test code to see how to connect with GPIB and query!!
#import visa
#rm = visa.ResourceManager()
#rm.list_resources()
#('ASRL1::INSTR','ASRL2::INSTR','GPIB0::19::INSTR')
#inst = rm.open_resource('GPIB0::19:INSTR')
#print(inst.query("*IDN?"))
#Ok Level 1 of testing complete - lets get started with the actual code now :) #
#Lets import all the packages we want to use in this program
#import numpy as np
#import matplotlib as plt
import visa, time, numpy as np
import csv
#Lets Initialize our resource manager to speak with the devices
rm=visa.ResourceManager()
#We want to do this just to check out spool results for reference later on.
rm.list_resources()
#Just checking if we are able to access both the instruments correctly.
agi4411b=rm.open_resource('GPIB0::18::INSTR')
agi8648a=rm.open_resource('GPIB0::19::INSTR')
print(agi4411b.query("*IDN?"))
print(agi8648a.query("*IDN?"))
frequency=raw_input('Please Input Frequency in MHZ\n')
power=raw_input('Please Input Power in dBm\n')
nooftrials=raw_input('Please Input Number of Trials for calibrating random error')
#Lets give the signal generator and spectrum analyzer their required basic settings.
agi8648a.write("OUTPUT 719;FREQ:CW %s MHZ" % frequency)
agi8648a.write("POW:AMPL %s DBM" % power)
agi4411b.write("SENS:FREQ:CENT %s MHZ" % frequency)
agi4411b.write("SENS:FREQ:SPAN 1 MHZ")
agi4411b.write("POW:ATT 0 DB")
agi8648a.write(":OUTP:STAT ON")
#Loop for selecting optimum frequency and amplitude: -
rbw_val = ['1 KHZ','10 KHZ','30 KHZ','100 KHZ','300 KHZ','1 MHZ']
j=0
while j<len(rbw_val):
    agi4411b.write(":SENSe:BWIDTH:RES ", rbw_val[j])# Write GPIB to select first RBW value and get the xval
    print("rbw set as ", rbw_val[j])
    agi4411b.write("INIT:IMM;*WAI")# to make a complete sweep before actually getting the data.
    agi4411b.write("CALC:MARK1:MAX")
    freqcstr=agi4411b.query("CALC:MARK1:X?")
    freqc=float(freqcstr)
    freqcchk=freqc/1e6
    freqchk=float(frequency)
    print(freqchk,freqcchk)
    if freqchk>freqcchk:
        j=j+1
        rbw_amp2=rbw_val[j]
        print("if loop j ",j)
        j= len(rbw_val)+1
        print(rbw_amp2)
        print("exiting calculation loop")
    elif freqchk<=freqcchk:
        j= j+1
        print(rbw_val[j])
        print("elif loop j ",j)
        rbw_amp1=rbw_val[j]
else:
    print("Optimized Spectrum Analyzer values calculated")
    print(rbw_amp1, rbw_amp2)
#Let us loop the loop for storing frequency and Amplitude.
i=0
notr=int(nooftrials)
result=[[nooftrials,frequency,power]]
while i < notr:
#for trial in nooftrials: using amp1 to calculate frequency for x value
    agi4411b.write(":SENSe:BWIDTH:RES ", rbw_val[0])#write GPIB to select RBW value for xval amp1
    agi4411b.write("INIT:IMM;*WAI")
    agi4411b.write("CALC:MARK1:MAX")
    xval=(agi4411b.query("CALC:MARK1:X?"))
    xval_fin=float(xval)/1e6
    agi4411b.write(":SENSe:BWIDTH:RES ", rbw_amp1)#write GPIB to select RBW value for xval amp1
    agi4411b.write("INIT:IMM;*WAI")
    agi4411b.write("CALC:MARK1:MAX")
    yval1=float(agi4411b.query("CALC:MARK1:Y?"))
    agi4411b.write(":SENSe:BWIDTH:RES ", rbw_amp2)#write GPIB to select RBW value for xval amp2
    agi4411b.write("INIT:IMM;*WAI")
    agi4411b.write("CALC:MARK1:MAX")
    yval2=float(agi4411b.query("CALC:MARK1:Y?"))
    print i,xval,yval1,yval2
    yval=np.average([yval1,yval2])
    result.append([i,xval_fin,yval])
    agi8648a.write(":OUTP:STAT OFF")
    time.sleep(1)
    agi8648a.write(":OUTP:STAT ON")
    time.sleep(2)
    i = i + 1
#    print i
#    print trial
else:
    print i, "Are the number of trials conducted successfully"
    print result
#save file to data for backup: -  with some sort of timestamp
    with open('output'+time.strftime('%Y%m%d-%H%M%S')+'.csv',"wb") as f:
        writer = csv.writer(f)
        writer.writerows(result)
#Close the instruments and plot the data.
agi8648a.write(":OUTP:STAT OFF")
agi4411b.close
agi8648a.close

max_stars_repo_name: kashchitta/ams-instrcontrol

max_issues_repo_name: kashchitta/ams-instrcontrol

max_forks_repo_name: kashchitta/ams-instrcontrol


text:
import sys

sys.path.insert(0, '..')
import json
import os
import torch
from torch import optim
import argparse
from torchvision import transforms
import numpy as np
from torch.utils.data import DataLoader
from common.helpers.video import Display
from common.helpers.train import *
from common.helpers.support import *
from slam.vision.visuals import FeatDisplay
from common.helpers.transfer_model import *
from lanes.dataset import LanesDataset
from lanes.model import LanesSegNet

if __name__ == '__main__':
    parser = argparse.ArgumentParser(
        description='AI over the self driving space.')
    parser.add_argument('--test',
                        type=bool,
                        default=False,
                        help='Trains Models')
    parser.add_argument('--batch',
                        type=int,
                        default=16,
                        help='Batch size of input')

    parser.add_argument(
        '--Model',
        type=int,
        default=0,
        help='Model #0: Lanes Net \n Model #1: Quantized Lanes Net')

    args = parser.parse_args()

    SEED = 98
    device = torch.device(
        'cuda') if torch.cuda.is_available() else torch.device('cpu')
    BATCH_SIZE = args.batch
    PIN_MEM = True
    NUM_WORKERS = 4
    IMG_SIZE = 512
    CLASSES = 4
    if not args.test:
        print("Train Mode.")
        EPOCHS = 100
    else:
        print("Test Mode.")
        EPOCHS = 1

    np.random.seed(SEED)
    torch.manual_seed(SEED)
    torch.backends.cudnn.benchmark = True

    # Preprocess Data
    preprocess = transforms.Compose([
        transforms.ToTensor(),
        transforms.Resize((IMG_SIZE, IMG_SIZE)),
    ])

    # Video Paths
    video_url = "../etc/videos/test.mp4"
    display = Display(video_url)

    # Loads paths
    json_file = open('paths.json')
    json_str = json_file.read()
    paths = json.loads(json_str)

    # Lane Segmentation
    loss_fn = torch.nn.MSELoss()
    model = LanesSegNet(CLASSES)
    optimizer = optim.Adam(model.parameters(), lr=1e-6)
    pt_display = FeatDisplay()
    pt_display.init()
    pt_display.show()
    # pt_display.test_vision()
    # pt_display.show_with_model()
    # trainset = LanesDataset(paths["img_path"], paths["lanes_mask_trainpath"],
    #                         preprocess, 6)
    # calibrated_set = LanesDataset(paths["img_path"],
    #                               paths["lanes_mask_trainpath"], preprocess,
    #                               36)
    # valset = LanesDataset(paths["img_val_path"], paths["lanes_mask_valpath"],
    #                       preprocess, 10)
    # calibrated_loader = DataLoader(calibrated_set,
    #                                BATCH_SIZE,
    #                                num_workers=NUM_WORKERS,
    #                                pin_memory=PIN_MEM,
    #                                shuffle=True)

    # train_loader = DataLoader(trainset,
    #                           BATCH_SIZE,
    #                           num_workers=NUM_WORKERS,
    #                           pin_memory=PIN_MEM,
    #                           shuffle=True)
    # val_loader = DataLoader(valset,
    #                         BATCH_SIZE,
    #                         num_workers=NUM_WORKERS,
    #                         pin_memory=PIN_MEM,
    #                         shuffle=False)

    # if not args.test:
    #     fit_model(model, optimizer, train_loader, loss_fn, device, EPOCHS,
    #               val_loader)
    # else:
    #     display.show(model, device)


text:
#!/usr/bin/env python3

# TODO from jump we need to check out_dir exists and we can create subdir for logs there [do we pass out_dir with args?]

# TODO how best to not clobber CSV output files

# TODO think through and implement where data in *general* gets saved and log file handling (rotation?)

# TODO can we write to Eric spreadsheet template directly?

# TODO improve initial prompt for go/no-go decision

# TODO all of the TSH commanding (set/get config is primary objective for now)

# TODO have TshAccelBuffer natively include some form of raw_data_from_socket and maybe selected signal proc. methods

# TODO implement new feature with ability to start somewhere besides +x and progress thru abbrev. sequence after that...
# TODO ...BUT this requires special initialization for position to *previous* rough home

# TODO add another buffer duration, so (1) for final output file (60sec), & new shorter one for (2) GSS angles/visits

# TODO fix the move_to code for robustness against not getting actual position to desired (motor on forcing)

# TODO add PGA and ELECT testing and output files

import os
import sys
import time
import datetime
import platform
import logging
import logging.config
import numpy as np
import multiprocessing

from tshcal.inputs import argparser
from tshcal.inputs import user_menu
from tshcal.commanding import tsh_commands
from tshcal.commanding import esp_commands
from tshcal.common import buffer
from tshcal.defaults import ROOT_DIR, DEFAULT_PORT
from tshcal.common.buffer import Tsh, raw_data_from_socket

def get_logger(log_file):
    """return logger object with both stream and file handlers, the latter using log_file"""

    # create logger
    logger = logging.getLogger('tshcal')
    logger.setLevel(logging.DEBUG)

    # create file handler which logs to DEBUG level
    fh = logging.FileHandler(log_file)
    fh.setLevel(logging.DEBUG)

    # # create console handler with a higher log level
    # ch = logging.StreamHandler()
    # ch.setLevel(logging.INFO)  # changed from ERROR to INFO

    # create formatter and add it to the handlers
    formatter = logging.Formatter('%(asctime)s - %(name)s - %(levelname)s - %(module)s : %(lineno)d - %(message)s')
    fh.setFormatter(formatter)
    # ch.setFormatter(formatter)

    # add the handlers to the logger
    logger.addHandler(fh)
    # logger.addHandler(ch)

    return logger

def get_inputs(module_logger):
    """return arguments parsed from command line and make log entries"""

    # parse arguments
    module_logger.debug('parsing input arguments')
    args = argparser.parse_inputs()
    module_logger.info(str(args).replace('Namespace', 'Inputs: '))

    # adjust log level based on verbosity input args
    if args.quiet:
        module_logger.warning('Now switching to quiet for logging, that is, log level = WARNING.')
        level = logging.WARNING
    elif args.verbose:
        level = logging.DEBUG
    else:
        level = logging.INFO

    # switch to same level for all handlers
    module_logger.setLevel(level)
    for handler in module_logger.handlers:
        handler.setLevel(level)

    return args

def wait_for_start_time(s, mod_logger):
    """delay until start time, s, has been surpassed; log via mod_logger"""

    # TODO make more user-friendly; e.g. prompt/update every 10sec with time remaining to start & a chance to abort

    # FIXME below lines are development/debug "if False" mechanism used so I am not waiting during code development
    # ----- WHEN READY GET RID OF if False PRETEXT AND UNINDENT THE OTHER LINES OF CODE TO DEPLOY ACTUAL DELAY -----
    if False:
        now = datetime.datetime.now()
        if s <= now:
            raise RuntimeError('not enough delay in start time, use more of a delay until start time')

        while now < s:
            time.sleep(10)  # wait 10 seconds
            now = datetime.datetime.now()

        mod_logger.info('The calibration start time, %s, has been reached.  Begin calibrating now.' % s)
    # ----- WHEN READY GET RID OF if False PRETEXT AND UNINDENT THE OTHER LINES OF CODE TO DEPLOY ACTUAL DELAY -----

    mod_logger.info('Faking that the calibration start time, %s, has been reached.  Begin calibrating now.' % s)

def show_tsh_buffer_summary(tsh, sec=3, logger=None):
    """this does not much beyond demonstrate instantiation of TshAccelBuffer with logging"""

    # FIXME significant parts of this function are for quick demo purposes only, so scrub and fix

    logger.info('Capturing quick data summary from tsh...')

    # create data buffer
    buff = buffer.TshAccelBuffer(tsh, sec, logger=logger)
    raw_data_from_socket(tsh.ip, buff, port=DEFAULT_PORT)  # this populates 2nd arg, buff

    logger.debug('Done capturing quick data summary from tsh.')

    # FIXME should we do something smart about units?
    s = 'QUICK SUMMARY...Mean: X = {:.3f}, Y = {:.3f}, Z = {:.3f} [units]'.format(*np.mean(buff.xyz, axis=0))
    s += ' <-' * 17

    # show/log summary
    logger.info(s)

def main():
    """return status/exit code that results from running tshcal main application"""

    # create logger
    log_file = os.path.join(ROOT_DIR, 'logging', 'tshcal.log')
    module_logger = get_logger(log_file)
    module_logger.info('- STARTING MAIN TSHCAL APP - - - - - - - - - - - - - - - - - - - - - -')
    module_logger.debug('log_file = %s' % log_file)

    # get input arguments
    args = get_inputs(module_logger)

    # prompt user to follow along with logging in new terminal
    prompt_str = 'Start new log term w/ "tail -f %s" to see logging, then back to this cmd term for prompts.' % log_file
    accept_str = 'User indicated log file is being monitored in another terminal, so continue.'
    ans = user_menu.prompt_yes_no(prompt_str, accept_str)
    # if user chose zero to quit early, then do a graceful exit
    if ans == 0:
        module_logger.info('bye')
        sys.exit(-1)

    # TODO design tsh class that gives robustness (with commanding to set/get sample rate, gain, units, etc.)
    # create tsh object FIXME << this is a dummy for now
    tsh = Tsh(args.sensor, args.rate, args.gain)

    # FIXME for now, just squawk about not having Tsh class to handle get/set commanding or querying state
    module_logger.info('SKIPPING TSH SET/GET SINCE NO GOOD Tsh CLASS YET.')
    if False:
        # FIXME this indented block under "if False" should eventually get unindented to actually be employed

        # FIXME verify the tsh object we created given desired input args matches what's physically hooked up

        # set tsh parameters
        tsh_state_desired = tsh_commands.set_tsh_state(args)

        # get tsh parameters
        tsh_state_actual = tsh_commands.get_tsh_state()

        # verify tsh parameters
        if tsh_state_actual == tsh_state_desired:
            module_logger.debug('The tsh actual state matches our desired state.')
        else:
            # TODO give more info here -- what exactly does not match?
            raise AssertionError('The tsh actual state does NOT match our desired state.')

    # if OS allows, then we can do quick tsh data check (must be linux)
    pform = platform.platform()
    if pform.lower().startswith('linux'):

        module_logger.info('Perform quick tsh data test since OS detected is linux.')

        # FIXME this next multiprocessing feature does not work on mac os
        # create buffer to capture 2 seconds of TSH data and show user a summary of what we got (do as process for timeout)
        buff_sec = 2
        p = multiprocessing.Process(target=show_tsh_buffer_summary,
                                    args=(tsh, ),
                                    kwargs={'sec': buff_sec, 'logger': module_logger})
        p.start()

        # wait for double buff_sec or until process finishes
        p.join(buff_sec * 2)

        # if thread still active, that is, active for too long, then we log it and quit out
        if p.is_alive():
            timeout_msg = 'A %d-sec tsh data buffer still not filled after %d sec...' % (buff_sec, buff_sec * 2)
            timeout_msg += 'Too long, something went wrong...kill it!'
            module_logger.info(timeout_msg)

            p.terminate()
            p.join()

            module_logger.info('Why did it take %d seconds or more to fill a %d-second TSH buffer?' %
                               (buff_sec * 2, buff_sec))

            print('See log for early exit due to tsh buffer issue.')

            sys.exit(-2)

    else:
        module_logger.info('Skip quick tsh data test since OS is not Linux.')

    # FIXME Do we need to do anything prep/config for ESP here? (e.g. GENERAL MODE SELECTION or STATUS FUNCTIONS...
    # FIXME ...maybe from Table 3.5.1 of ESP301 User Guide or possibly something else)?  Will may have answered this?

    # FIXME before any motion begins, our safe trajectories ASSUME starting in home position; prompt user to verify home

    # FIXME Before delay timer, prompt user with enough info here to enable a go/no-go answer...
    # FIXME ... e.g. hit [Enter] if prompted-profile looks good OR give another keypress with chance to retry or...?
    # FIXME ** This should include whether or not debug/prompt-along-way is going to happen or not [auto- or semi-auto?]

    # offer enough prompt info for user to make go/no-go decision...but for now, just this crude prompt
    prompt_str = 'ESP rig at absolute home in all 3 axes (and other stuff initialized)?'
    accept_str = 'User indicated that ESP rig, TSH, cable clearance, config and all is a "go".'
    ans_go = user_menu.prompt_yes_no(prompt_str, accept_str)
    # if user chose zero to quit early, then do a graceful exit
    if ans_go == 0:
        module_logger.info('User aborted at go/no-go decision prompt.')
        sys.exit(-3)

    # FIXME our flow chart shows 2 delays, but no smarts here yet to verify enough time for TSH temperature settling

    # delay until start time to begin calibration
    wait_for_start_time(args.start, module_logger)

    # run calibration routine
    esp_commands.run_cal(tsh, args.outdir, plot=args.plot, debug=args.debug)

    # FIXME are there any commands we need to send to TSH at this point after running calibration?

    module_logger.info('- - - Calibration Complete - - - - - - - - - - - - - - - - - - - - - -')

    return 0  # return zero for success (typical Linux command line behavior)

if __name__ == '__main__':

    sys.exit(main())

max_stars_repo_name: kenhro/tshcal

max_issues_repo_name: kenhro/tshcal

max_forks_repo_name: kenhro/tshcal


text:
#===============================================================================
# acoular example1.py without frequency beamformers
#===============================================================================

from os import path
import acoular
from pylab import figure, imshow, colorbar, show, subplot,title
from numpy import zeros

micgeofile = path.join(path.split(acoular.__file__)[0],'xml','array_56.xml')
datafile = 'example_data.h5'
calibfile = 'example_calib.xml'
cfreq = 4000

mg = acoular.MicGeom( from_file=micgeofile )
ts = acoular.MaskedTimeSamples( name=datafile )

ts.start = 0
ts.stop = 16000
invalid = [1,7]
ts.invalid_channels = invalid
mg.invalid_channels = invalid
rg = acoular.RectGrid( x_min=-0.6, x_max=0.0, y_min=-0.3, y_max=0.3, z=0.68, \
increment=0.05 )
bt = acoular.BeamformerTime(source = ts, grid = rg, mpos = mg, c=346.04)
ft = acoular.FiltFiltOctave(source=bt, band=cfreq)
pt = acoular.TimePower(source=ft)
avgt = acoular.TimeAverage(source=pt, naverage = 1024)
cacht = acoular.TimeCache( source = avgt) # cache to prevent recalculation

#===============================================================================
# delay and sum beamformer in time domain with autocorrelation removal
# processing chain: zero-phase filtering, beamforming+power, average
#===============================================================================

fi = acoular.FiltFiltOctave(source=ts, band=cfreq)
bts = acoular.BeamformerTimeSq(source = fi,grid=rg, mpos=mg, r_diag=True,c=346.04)
avgts = acoular.TimeAverage(source=bts, naverage = 1024)
cachts = acoular.TimeCache( source = avgts) # cache to prevent recalculation

#===============================================================================
# plot result maps for different beamformers in time domain
#===============================================================================
i1 = 1
i2 = 2 # no of figure
for block in (cacht, cachts):
    # first, plot time-dependent result (block-wise)
    figure(i2)
    i2 += 1
    res = zeros(rg.size) # init accumulator for average
    i3 = 1 # no of subplot
    for block_res in block.result(1):  #one single block
        subplot(4,4,i3)
        i3 += 1
        res += block_res[0] # average accum.
        map = block_res[0].reshape(rg.shape)
        mx = acoular.L_p(map.max())
        imshow(acoular.L_p(map.T), vmax=mx, vmin=mx-15,
               interpolation='nearest', extent=rg.extend())
        title('%i' % ((i3-1)*1024))
    res /= i3-1 # average
    # second, plot overall result (average over all blocks)
    figure(1)
    subplot(3,4,i1)
    i1 += 1
    map = res.reshape(rg.shape)
    mx = acoular.L_p(map.max())
    imshow(acoular.L_p(map.T), vmax=mx, vmin=mx-15,
           interpolation='nearest', extent=rg.extend())
    colorbar()
    title(('BeamformerTime','BeamformerTimeSq')[i2-3])
show()


text:
import numpy as np
import utils
from torch import nn

def validation_binary(model: nn.Module, criterion, valid_loader):
    model.eval()
    losses = []

    jaccard = []
    dice = []

    for inputs, targets in valid_loader:
        inputs = utils.variable(inputs, volatile=True)
        targets = utils.variable(targets)
        outputs = model(inputs)
        loss = criterion(outputs, targets)
        losses.append(loss.item())
        outputs = (outputs > 0).float()
        jaccard += get_jaccard(targets, outputs)
        dice += get_dice(targets, outputs)

    valid_loss = np.mean(losses)  # type: float

    valid_jaccard = np.mean(jaccard).astype(np.float64)
    valid_dice = np.mean(dice).astype(np.float64)

    print('Valid loss: {:.5f}, jaccard: {:.5f}, dice: {:.5f}'.format(valid_loss, valid_jaccard, valid_dice))
    metrics = {'valid_loss': valid_loss, 'jaccard_loss': valid_jaccard, 'dice_loss': valid_dice}
    return metrics

def multi_validation_binary(model: nn.Module, criterion, valid_loader):
    model.eval()
    losses = []
    slosses = []
    closses = []

    jaccard = []
    dice = []

    for inputs, targets, ctargets in valid_loader:
        inputs = utils.variable(inputs, volatile=True)
        targets = utils.variable(targets)
        ctargets = utils.variable(ctargets)
        outputs, coutputs = model(inputs)
        sloss, closs = criterion(outputs, targets, coutputs, ctargets)
        loss = sloss + closs
        losses.append(loss.item())
        slosses.append(sloss.item())
        closses.append(closs.item())
        outputs = (outputs > 0).float()
        jaccard += get_jaccard(targets, outputs)
        dice += get_dice(targets, outputs)

    valid_loss = np.mean(losses)  # type: float
    svalid_loss = np.mean(slosses)
    cvalid_loss = np.mean(closses)

    valid_jaccard = np.mean(jaccard).astype(np.float64)
    valid_dice = np.mean(dice).astype(np.float64)

    print('Valid loss: {:.5f}, jaccard: {:.5f}, dice: {:.5f}'.format(valid_loss, valid_jaccard, valid_dice))
    metrics = {'valid_loss': valid_loss, 's_valid_loss': svalid_loss, 'c_valid_loss': cvalid_loss, 'jaccard_loss': valid_jaccard, 'dice_loss': valid_dice}
    return metrics

def get_jaccard(y_true, y_pred):
    epsilon = 1e-15
    intersection = (y_pred * y_true).sum(dim=-2).sum(dim=-1)
    union = y_true.sum(dim=-2).sum(dim=-1) + y_pred.sum(dim=-2).sum(dim=-1)

    return list((intersection / (union + epsilon - intersection)).data.cpu().numpy())

def get_dice(y_true, y_pred):
    epsilon = 1e-15
    intersection = (y_pred * y_true).sum(dim=-2).sum(dim=-1)
    union = y_true.sum(dim=-2).sum(dim=-1) + y_pred.sum(dim=-2).sum(dim=-1)
    return list((2 * intersection / (union + epsilon)).data.cpu().numpy())

max_stars_repo_name: crossknight/angiodysplasia-segmentation

max_issues_repo_name: crossknight/angiodysplasia-segmentation

max_forks_repo_name: crossknight/angiodysplasia-segmentation


text:
import numpy as np
import SimpleITK as sitk
from pathlib import Path

#normalization
def whitening(image):
    """Whitening. Normalises image to zero mean and unit variance."""

    image = image.astype(np.float32)

    mean = np.mean(image)
    std = np.std(image)

    if std > 0:
        ret = (image - mean) / std
    else:
        ret = image * 0.
    return ret

def normalise_zero_one(image):
    """Image normalisation. Normalises image to fit [0, 1] range."""

    image = image.astype(np.float32)

    minimum = np.min(image)
    maximum = np.max(image)

    if maximum > minimum:
        ret = (image - minimum) / (maximum - minimum)
    else:
        ret = image * 0.
    return ret

def normalise_one_one(image):
    """Image normalisation. Normalises image to fit [-1, 1] range."""

    ret = normalise_zero_one(image)
    ret *= 2.
    ret -= 1.
    return ret


text:
import json
import random
import numpy as np

# available moves
moves = ['rock', 'paper', 'scissors', 'lizard', 'spock']
# adjacency matrix for game
game = [[ 0, -1,  1,  1, -1],
        [ 1,  0, -1, -1,  1],
        [-1,  1,  0,  1, -1],
        [-1,  1, -1,  0,  1],
        [ 1, -1,  1, -1,  0]]

# model data
# this is what should be saved to CosmosDb
model = {
    'hist': 5,
    'play': [],
    'pred': []
}

# this can either be retrieved from
# CosmosDb or 0-initialized as belo
def init():
    for i in range(model['hist']):
        if i == 0:
            model['pred'].append({
                'count': 0,
                'probs': { }
            })
        else:
            model['pred'].append({ })

def predict():
    history = model['play']
    pred = model['pred']
    guesses = np.zeros((model['hist'], len(moves)))

    for i in range(model['hist']):
        record = pred[i]

        # guessing based on user distribution
        if i == 0:
            for k in record['probs'].keys():
                idx = moves.index(k)
                guesses[i][idx] = record['probs'][k] / record['count']
        # guess on conditional distribution
        else:
            m = history[-i:]
            # not enough history
            if len(m) < i:
                break

            item = '_'.join(m)
            if item in record:
                for k in record[item]['probs'].keys():
                    idx = moves.index(k)
                    guesses[i][idx] = record[item]['probs'][k] / record[item]['count']

    # nothing to go on so just guess
    if np.sum(guesses) == 0:
        return random.choice(moves)
    else:
        # get the highest probability slots
        probabilities = np.max(guesses, axis=0)
        max_p = np.max(probabilities)
        # may have more than one max_p
        their_moves = [idx for idx, val in enumerate(probabilities) if val == max_p]

        # fill out their choices
        choices = []
        for guess in their_moves:
            # find out moves that beats projected move
            choices.extend([moves[idx] for idx, val in enumerate(game[guess]) if val == -1])

        # choose one randomly
        return random.choice(choices)

def update(move):
    history = model['play']
    predictions = model['pred']

    # add to history
    history.append(move)

    # update model
    for i in range(model['hist']):
        # last moves i moves
        m = history[-(i+1):]

        # if history is insufficient
        if len(m) == i:
            break

        # get conditional
        key = m[:i]
        # get vals
        val = m[-1]

        pred = predictions[i]

        # update predictions n_grams
        if i == 0:
            predictions[i]['count'] += 1
            if val in pred['probs']:
                pred['probs'][val] += 1
            else:
                pred['probs'][val] = 1
        else:
            skey = '_'.join(key)
            if skey in pred:
                pred[skey]['count'] += 1
                if val in pred[skey]['probs']:
                    pred[skey]['probs'][val] += 1
                else:
                    pred[skey]['probs'][val] = 1
            else:
                pred[skey] = {
                    'count': 1,
                    'probs': { val : 1 }
                }

def winner(player1, player2):
    idx1, idx2 = moves.index(player1), moves.index(player2)
    w = game[idx1][idx2]
    return w

if __name__ == "__main__":
    # init model
    init()

    plays = ['rock', 'paper', 'scissors', 'rock', 'paper', 'scissors', 'lizard', 'rock', 'spock', 'rock', 'paper']

    # fill model (fast forward through plays)
    for item in plays:
        update(item)

    print(f'Model:\n{json.dumps(model, sort_keys=True,indent=3)}\n')

    for i in range(len(plays)):
        print(f'\nRound {i+1}')
        # computer guess
        computer = predict()
        # human move
        human = plays[i]

        print(f'You: {human} vs. Computer: {computer}')

        # who won??
        w = winner(human, computer)
        print(f'{["computer wins", "a tie", "human wins"][w+1]}')

        # update model
        update(human)

max_stars_repo_name: sethjuarez/gamemoves

max_issues_repo_name: sethjuarez/gamemoves

max_forks_repo_name: sethjuarez/gamemoves


text:
import numpy as np
import cv2

# VIDEO = "./data/owata_ta.mkv"
# VIDEO = "./data/a.flv"
# VIDEO = "./data/unsaga.mp4"
VIDEO = 0
cap = cv2.VideoCapture(VIDEO)

kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(3,3))
fgbg = cv2.bgsegm.createBackgroundSubtractorLSBP(
    nSamples=20,
    LSBPRadius=16,
    Tlower=2.0,
    Tupper=32.0,
    Tinc=1.0,
    Tdec=0.05,
    Rscale=10.0,
    Rincdec=0.005,
    LSBPthreshold=8,
)

while(1):
    ret, frame = cap.read()
    if not ret:
        break
    fgmask = fgbg.apply(frame)
    # fgmask = cv2.morphologyEx(fgmask, cv2.MORPH_OPEN, kernel)

    cv2.imshow('frame',fgmask)
    k = cv2.waitKey(30) & 0xff
    if k == 27:
        break

cap.release()
cv2.destroyAllWindows()


text:
# -*- coding: utf-8 -*-
from aiida_quantumespresso.calculations.pw import PwCalculation
from aiida_quantumespresso.parsers.basic_raw_parser_pw import (
    parse_raw_output, QEOutputParsingError)
from aiida.orm.data.parameter import ParameterData
from aiida.orm.data.folder import FolderData
from aiida.parsers.parser import Parser  # , ParserParamManager
from aiida_quantumespresso.parsers import convert_qe2aiida_structure
from aiida.common.datastructures import calc_states
from aiida.common.exceptions import UniquenessError
from aiida.orm.data.array import ArrayData
from aiida.orm.data.array.kpoints import KpointsData

# TODO: I don't like the generic class always returning a name for the link to the output structure

class BasicpwParser(Parser):
    """
    This class is the implementation of the Parser class for PWscf.
    """

    _setting_key = 'parser_options'

    def __init__(self, calc):
        """
        Initialize the instance of PwParser
        """
        # check for valid input
        if not isinstance(calc, PwCalculation):
            raise QEOutputParsingError("Input calc must be a PwCalculation")

        super(BasicpwParser, self).__init__(calc)

    def parse_with_retrieved(self, retrieved):
        """
        Receives in input a dictionary of retrieved nodes.
        Does all the logic here.
        """
        from aiida.common.exceptions import InvalidOperation
        import os
        import glob

        successful = True

        # retrieve the input parameter
        calc_input = self._calc.inp.parameters

        # look for eventual flags of the parser
        try:
            parser_opts = self._calc.inp.settings.get_dict()[self.get_parser_settings_key()]
        except (AttributeError, KeyError):
            parser_opts = {}

        # load the input dictionary
        # TODO: pass this input_dict to the parser. It might need it.
        input_dict = self._calc.inp.parameters.get_dict()

        # Check that the retrieved folder is there
        try:
            out_folder = retrieved[self._calc._get_linkname_retrieved()]
        except KeyError:
            self.logger.error("No retrieved folder found")
            return False, ()

        # check what is inside the folder
        list_of_files = out_folder.get_folder_list()
        # at least the stdout should exist
        if not self._calc._OUTPUT_FILE_NAME in list_of_files:
            self.logger.error("Standard output not found")
            successful = False
            return successful, ()
        # if there is something more, I note it down, so to call the raw parser
        # with the right options
        # look for xml
        has_xml = False
        if self._calc._DATAFILE_XML_BASENAME in list_of_files:
            has_xml = True
        # look for bands
        has_bands = False
        if glob.glob(os.path.join(out_folder.get_abs_path('.'),
                                  'K*[0-9]')):
            # Note: assuming format of kpoints subfolder is K*[0-9]
            has_bands = True
            # TODO: maybe it can be more general than bands only?
        out_file = os.path.join(out_folder.get_abs_path('.'),
                                self._calc._OUTPUT_FILE_NAME)
        xml_file = os.path.join(out_folder.get_abs_path('.'),
                                self._calc._DATAFILE_XML_BASENAME)
        dir_with_bands = out_folder.get_abs_path('.')

        # call the raw parsing function
        parsing_args = [out_file, input_dict, parser_opts]
        if has_xml:
            parsing_args.append(xml_file)
        if has_bands:
            if not has_xml:
                self.logger.warning("Cannot parse bands if xml file not "
                                    "found")
            else:
                parsing_args.append(dir_with_bands)

        out_dict, trajectory_data, structure_data, raw_successful = parse_raw_output(*parsing_args)

        # if calculation was not considered failed already, use the new value
        successful = raw_successful if successful else successful

        new_nodes_list = []

        # I eventually save the new structure. structure_data is unnecessary after this
        in_struc = self._calc.get_inputs_dict()['structure']
        type_calc = input_dict['CONTROL']['calculation']
        struc = in_struc
        if type_calc in ['relax', 'vc-relax', 'md', 'vc-md']:
            if 'cell' in structure_data.keys():
                struc = convert_qe2aiida_structure(structure_data, input_structure=in_struc)
                new_nodes_list.append((self.get_linkname_outstructure(), struc))

        k_points_list = trajectory_data.pop('k_points', None)
        k_points_weights_list = trajectory_data.pop('k_points_weights', None)

        if k_points_list is not None:

            # build the kpoints object
            if out_dict['k_points_units'] not in ['2 pi / Angstrom']:
                raise QEOutputParsingError('Error in kpoints units (should be cartesian)')
            # converting bands into a BandsData object (including the kpoints)

            kpoints_from_output = KpointsData()
            kpoints_from_output.set_cell_from_structure(struc)
            kpoints_from_output.set_kpoints(k_points_list, cartesian=True,
                                            weights=k_points_weights_list)
            kpoints_from_input = self._calc.inp.kpoints
            try:
                kpoints_from_input.get_kpoints()
            except AttributeError:
                new_nodes_list += [(self.get_linkname_out_kpoints(), kpoints_from_output)]

        # convert the dictionary into an AiiDA object
        output_params = ParameterData(dict=out_dict)
        # return it to the execmanager
        new_nodes_list.append((self.get_linkname_outparams(), output_params))

        if trajectory_data:
            import numpy
            from aiida.orm.data.array.trajectory import TrajectoryData
            from aiida.orm.data.array import ArrayData

            try:
                positions = numpy.array(trajectory_data.pop('atomic_positions_relax'))
                try:
                    cells = numpy.array(trajectory_data.pop('lattice_vectors_relax'))
                    # if KeyError, the MD was at fixed cell
                except KeyError:
                    cells = numpy.array([in_struc.cell] * len(positions))

                symbols = numpy.array([str(i.kind_name) for i in in_struc.sites])
                stepids = numpy.arange(len(positions))  # a growing integer per step
                # I will insert time parsing when they fix their issues about time
                # printing (logic is broken if restart is on)

                traj = TrajectoryData()
                traj.set_trajectory(stepids=stepids,
                                    cells=cells,
                                    symbols=symbols,
                                    positions=positions,
                )
                for x in trajectory_data.iteritems():
                    traj.set_array(x[0], numpy.array(x[1]))
                # return it to the execmanager
                new_nodes_list.append((self.get_linkname_outtrajectory(), traj))

            except KeyError:  # forces in scf calculation (when outputed)
                arraydata = ArrayData()
                for x in trajectory_data.iteritems():
                    arraydata.set_array(x[0], numpy.array(x[1]))
                # return it to the execmanager
                new_nodes_list.append((self.get_linkname_outarray(), arraydata))

        return successful, new_nodes_list

    def get_parser_settings_key(self):
        """
        Return the name of the key to be used in the calculation settings, that
        contains the dictionary with the parser_options
        """
        return 'parser_options'

    def get_linkname_outstructure(self):
        """
        Returns the name of the link to the output_structure
        Node exists if positions or cell changed.
        """
        return 'output_structure'

    def get_linkname_outtrajectory(self):
        """
        Returns the name of the link to the output_trajectory.
        Node exists in case of calculation='md', 'vc-md', 'relax', 'vc-relax'
        """
        return 'output_trajectory'

    def get_linkname_outarray(self):
        """
        Returns the name of the link to the output_array
        Node may exist in case of calculation='scf'
        """
        return 'output_array'

    def get_linkname_out_kpoints(self):
        """
        Returns the name of the link to the output_kpoints
        Node exists if cell has changed and no bands are stored.
        """
        return 'output_kpoints'

max_stars_repo_name: unkcpz/aiida-quantumespresso

max_issues_repo_name: unkcpz/aiida-quantumespresso

max_forks_repo_name: unkcpz/aiida-quantumespresso


text:
import os
import numpy as np

'''
Used to separate the road images alone from final_data.npy
'''

final_data_dir = "C:\\Users\\gauta\\python projects\\self-driving-car\\data\\final_data.npy"
org_data = np.load(final_data_dir, allow_pickle=True)

road_data = np.array([i for i, _, _ in org_data])
print(road_data.shape)

if not os.path.exists('data'):
    os.makedirs('data')

np.save('data\\road_data.npy', road_data)
print('Road data saved!')


text:
# OBS
# Imagem celular original do vírus recentemente descoberto SARS-CoV-2,
# popularmente chamado de COVID-19 ou Coronavirus.

import cv2 as cv
import numpy as np
import matplotlib.pyplot as plt

!wget "https://raw.githubusercontent.com/PedroHaupenthal/Image-Processing/master/watershed/covid_19.jpg" -O "covid_19.jpg"

img1 = cv.imread("covid_19.jpg")
img1 = cv.cvtColor(img1, cv.COLOR_BGR2RGB)

img2 = cv.cvtColor(img1, cv.COLOR_RGB2GRAY)
img2 = cv.bitwise_not(img2)
ret, img2 = cv.threshold(img2,
                         0,
                         255,
                         cv.THRESH_BINARY_INV + cv.THRESH_OTSU)

kernel = np.ones((3,3), np.uint8)

img3 = cv.morphologyEx(img2,
                       cv.MORPH_CLOSE,
                       kernel,
                       iterations = 2)
img4 = cv.dilate(img3,
                 kernel,
                 iterations = 5)

img5 = cv.distanceTransform(img3,
                            cv.DIST_L2,
                            5)
ret,img6 = cv.threshold(img5,
                        0.65 * img5.max(),
                        255,
                        0)

img6 = np.uint8(img6)
img7 = cv.subtract(img4, img6)

ret, count = cv.connectedComponents(img6)
count = count + 1

count[img7 == 255] = 0
img8 = cv.watershed(img1, count)
img1[count == -1] = [255, 0, 0]

plt.figure(figsize=(30,30))
plt.subplot(121), plt.imshow(img1), plt.title("ORIGINAL"), plt.axis("off")
plt.subplot(122), plt.imshow(img8, cmap='jet'), plt.title("RESULTADO"), plt.axis("off")
plt.show()

max_stars_repo_name: PedroHaupenthal/image-processing

max_issues_repo_name: PedroHaupenthal/image-processing

max_forks_repo_name: PedroHaupenthal/image-processing