text:
# Copyright (c) 2016, Konstantinos Kamnitsas
# All rights reserved.
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the BSD license. See the accompanying LICENSE file
# or read the terms at https://opensource.org/licenses/BSD-3-Clause.

from __future__ import absolute_import, print_function, division

import os
import sys
import time
import numpy as np
import math
import random
import traceback
import multiprocessing
import signal
import collections

from deepmedic.dataManagement.io import loadVolume
from deepmedic.dataManagement.preprocessing import calculateTheZeroIntensityOf3dImage, padCnnInputs
from deepmedic.neuralnet.pathwayTypes import PathwayTypes as pt
from deepmedic.dataManagement.augmentation import augment_patch

# Order of calls:
# getSampledDataAndLabelsForSubepoch
#    get_random_subjects_to_train_subep
#    get_num_samples_per_cat_per_subj
#    load_subj_and_get_samples
#        load_imgs_of_subject
#        sample_coords_of_segments
#        extractSegmentGivenSliceCoords
#            getImagePartFromSubsampledImageForTraining
#    shuffleSegmentsOfThisSubepoch

# Main sampling process during training. Executed in parallel while training on a batch on GPU.
# Called from training.do_training()
def getSampledDataAndLabelsForSubepoch( log,
                                        train_or_val,
                                        num_parallel_proc,
                                        run_input_checks,
                                        cnn3d,
                                        max_n_cases_per_subep,
                                        n_samples_per_subep,
                                        sampling_type,

                                        paths_per_chan_per_subj,
                                        paths_to_lbls_per_subj,

                                        provided_mask,
                                        paths_to_masks_per_subj,

                                        provided_weightmaps_per_sampl_cat,
                                        paths_to_wmaps_per_sampl_cat_per_subj,

                                        pad_input_imgs,
                                        augm_params
                                        ):
    # Returns: channelsOfSegmentsForSubepochArrayPerPathway - List of arrays, one per pathway. Each array [N_samples, Channs, R,C,Z]
    #        lblsForPredictedPartOfSegmentsForSubepochArray - Array of shape: [N_samples, R_out, C_out, Z_out)

    id_str = "[SAMPLER-TR|PID:"+str(os.getpid())+"]" if train_or_val == "train" else "[SAMPLER-VAL|PID:"+str(os.getpid())+"]"
    start_time_sampling = time.time()
    training_or_validation_str = "Training" if train_or_val == "train" else "Validation"

    log.print3(id_str+" :=:=:=:=:=:=: Starting to sample for next [" + training_or_validation_str + "]... :=:=:=:=:=:=:")

    total_number_of_subjects = len(paths_per_chan_per_subj)
    inds_of_subjects_for_subep = get_random_subjects_to_train_subep( total_number_of_subjects = total_number_of_subjects,
                                                                     max_subjects_on_gpu_for_subepoch = max_n_cases_per_subep,
                                                                     get_max_subjects_for_gpu_even_if_total_less = False )
    log.print3(id_str+" Out of [" + str(total_number_of_subjects) + "] subjects given for [" + training_or_validation_str + "], "+
               "we will sample from maximum [" + str(max_n_cases_per_subep) + "] per subepoch.")
    log.print3(id_str+" Shuffled indices of subjects that were randomly chosen: "+str(inds_of_subjects_for_subep))

    # List, with [numberOfPathwaysThatTakeInput] sublists. Each sublist is list of [partImagesLoadedPerSubepoch] arrays [channels, R,C,Z].
    channelsOfSegmentsForSubepochListPerPathway = [ [] for i in range(cnn3d.getNumPathwaysThatRequireInput()) ]
    lblsForPredictedPartOfSegmentsForSubepochList = [] # Labels only for the central/predicted part of segments.
    n_subjects_for_subep = len(inds_of_subjects_for_subep) #Can be different than max_n_cases_per_subep, cause of available images number.

    # This is to separate each sampling category (fore/background, uniform, full-image, weighted-classes)
    n_samples_per_cat_per_subj = get_num_samples_per_cat_per_subj(n_samples_per_subep,
                                                                  sampling_type.getPercentOfSamplesPerCategoryToSample(),
                                                                  n_subjects_for_subep )

    args_sampling_job = [log,
                        train_or_val,
                        run_input_checks,
                        cnn3d,
                        sampling_type,
                        paths_per_chan_per_subj,
                        paths_to_lbls_per_subj,
                        provided_mask,
                        paths_to_masks_per_subj,
                        provided_weightmaps_per_sampl_cat,
                        paths_to_wmaps_per_sampl_cat_per_subj,
                        # Pre-processing:
                        pad_input_imgs,
                        augm_params,

                        n_subjects_for_subep,
                        inds_of_subjects_for_subep,
                        n_samples_per_cat_per_subj ]

    log.print3(id_str+" Will sample from [" + str(n_subjects_for_subep) + "] subjects for next " + training_or_validation_str + "...")

    jobs_inds_to_do = list(range(n_subjects_for_subep)) # One job per subject.

    if num_parallel_proc <= 0: # Sequentially
        for job_i in jobs_inds_to_do :
            (channelsOfSegmentsFromThisSubjPerPathway,
            lblsOfPredictedPartOfSegmentsFromThisSubj) = load_subj_and_get_samples( *( [job_i]+args_sampling_job ) )
            for pathway_i in range(cnn3d.getNumPathwaysThatRequireInput()) :
                channelsOfSegmentsForSubepochListPerPathway[pathway_i] += channelsOfSegmentsFromThisSubjPerPathway[pathway_i] # concat does not copy.
            lblsForPredictedPartOfSegmentsForSubepochList += lblsOfPredictedPartOfSegmentsFromThisSubj # concat does not copy.

    else: # Parallelize sampling from each subject
        while len(jobs_inds_to_do) > 0: # While jobs remain.
            jobs = collections.OrderedDict()

            log.print3(id_str+" ******* Spawning children processes to sample from [" + str(len(jobs_inds_to_do)) + "] subjects*******")
            log.print3(id_str+" MULTIPROC: Number of CPUs detected: " + str(multiprocessing.cpu_count()) + ". Requested to use max: [" + str(num_parallel_proc) + "]")
            num_workers = min(num_parallel_proc, multiprocessing.cpu_count())
            log.print3(id_str+" MULTIPROC: Spawning [" + str(num_workers) + "] processes to load data and sample.")
            worker_pool = multiprocessing.Pool(processes=num_workers, initializer=init_sampling_proc)

            try: # Stacktrace in multiprocessing: https://jichu4n.com/posts/python-multiprocessing-and-exceptions/
                for job_i in jobs_inds_to_do: # submit jobs
                    jobs[job_i] = worker_pool.apply_async( load_subj_and_get_samples, ( [job_i]+args_sampling_job ) )

                for job_i in list(jobs_inds_to_do): # copy with list(...), so that this loops normally even if something is removed from list.
                    try:
                        (channelsOfSegmentsFromThisSubjPerPathway,
                        lblsOfPredictedPartOfSegmentsFromThisSubj) = jobs[job_i].get(timeout=10) # timeout in case process for some reason never started (happens in py3)
                        for pathway_i in range(cnn3d.getNumPathwaysThatRequireInput()) :
                            channelsOfSegmentsForSubepochListPerPathway[pathway_i] += channelsOfSegmentsFromThisSubjPerPathway[pathway_i] # concat does not copy.
                        lblsForPredictedPartOfSegmentsForSubepochList += lblsOfPredictedPartOfSegmentsFromThisSubj # concat does not copy.
                        jobs_inds_to_do.remove(job_i)
                    except multiprocessing.TimeoutError as e:
                        log.print3(id_str+"\n\n WARN: MULTIPROC: Caught TimeoutError when getting results of job [" + str(job_i) + "].\n" +
                                   " WARN: MULTIPROC: Will resubmit job [" + str(job_i) + "].\n")
                        if num_workers == 1:
                            break # If this 1 worker got stuck, every job will wait timeout. Slow. Recreate pool.

            except (Exception, KeyboardInterrupt) as e:
                log.print3(id_str+"\n\n ERROR: Caught exception in getSampledDataAndLabelsForSubepoch(): "+str(e)+"\n")
                log.print3( traceback.format_exc() )
                worker_pool.terminate()
                worker_pool.join() # Will wait. A KeybInt will kill this (py3)
                raise e
            except: # Catches everything, even a sys.exit(1) exception.
                log.print3(id_str+"\n\n ERROR: Unexpected error in getSampledDataAndLabelsForSubepoch(). System info: ", sys.exc_info()[0])
                worker_pool.terminate()
                worker_pool.join()
                raise Exception("Unexpected error.")
            else: # Nothing went wrong
                worker_pool.terminate() # # Needed in case any processes are hunging. worker_pool.close() does not solve this.
                worker_pool.join()

    # Got all samples for subepoch. Now shuffle them, together segments and their labels.
    [channelsOfSegmentsForSubepochListPerPathway,
    lblsForPredictedPartOfSegmentsForSubepochList ] = shuffleSegmentsOfThisSubepoch(channelsOfSegmentsForSubepochListPerPathway,
                                                                                lblsForPredictedPartOfSegmentsForSubepochList )
    end_time_sampling = time.time()
    log.print3(id_str+" TIMING: Sampling for next [" + training_or_validation_str + "] lasted: {0:.1f}".format(end_time_sampling-start_time_sampling)+" secs.")

    log.print3(id_str+" :=:=:=:=:=:=: Finished sampling for next [" + training_or_validation_str + "] :=:=:=:=:=:=:")

    channelsOfSegmentsForSubepochArrayPerPathway = [ np.asarray(imPartsForPathwayi, dtype="float32") for imPartsForPathwayi in channelsOfSegmentsForSubepochListPerPathway ]
    lblsForPredictedPartOfSegmentsForSubepochArray = np.asarray(lblsForPredictedPartOfSegmentsForSubepochList, dtype="int32") # Could be int16 to save RAM?

    return [channelsOfSegmentsForSubepochArrayPerPathway,
            lblsForPredictedPartOfSegmentsForSubepochArray]

def init_sampling_proc():
    # This will make child-processes ignore the KeyboardInterupt (sigInt). Parent will handle it.
    # See: http://stackoverflow.com/questions/11312525/catch-ctrlc-sigint-and-exit-multiprocesses-gracefully-in-python/35134329#35134329
    signal.signal(signal.SIGINT, signal.SIG_IGN)

def get_random_subjects_to_train_subep( total_number_of_subjects,
                                        max_subjects_on_gpu_for_subepoch,
                                        get_max_subjects_for_gpu_even_if_total_less=False ):
    # Returns: list of indices
    subjects_indices = list(range(total_number_of_subjects)) #list() for python3 compatibility, as range cannot get assignment in shuffle()
    random_order_chosen_subjects=[]
    random.shuffle(subjects_indices) #does it in place. Now they are shuffled

    if max_subjects_on_gpu_for_subepoch>=total_number_of_subjects:
        random_order_chosen_subjects += subjects_indices

        if get_max_subjects_for_gpu_even_if_total_less : #This is if I want to have a certain amount on GPU, even if total subjects are less.
            while (len(random_order_chosen_subjects)<max_subjects_on_gpu_for_subepoch):
                random.shuffle(subjects_indices)
                number_of_extra_subjects_to_get_to_fill_gpu = min(max_subjects_on_gpu_for_subepoch - len(random_order_chosen_subjects), total_number_of_subjects)
                random_order_chosen_subjects += (subjects_indices[:number_of_extra_subjects_to_get_to_fill_gpu])
            assert len(random_order_chosen_subjects) != max_subjects_on_gpu_for_subepoch
    else:
        random_order_chosen_subjects += subjects_indices[:max_subjects_on_gpu_for_subepoch]

    return random_order_chosen_subjects

def get_num_samples_per_cat_per_subj(n_samples_per_subep,
                                     perc_samples_per_cat,
                                     n_subjects_for_subep ) :
    # perc_samples_per_cat: list with a percentage for each type of category to sample.

    # First, distribute n_samples_per_subep to the categories:
    n_sampl_categs = len(perc_samples_per_cat)
    n_samples_per_cat = np.zeros( n_sampl_categs, dtype="int32" )
    for cat_i in range(n_sampl_categs) :
        n_samples_per_cat[cat_i] += int(n_samples_per_subep * perc_samples_per_cat[cat_i])
    # Distribute samples that were left if perc dont exactly add to 1.
    n_undistributed_samples = n_samples_per_subep - np.sum(n_samples_per_cat)
    cats_to_distribute_samples = np.random.choice(n_sampl_categs, size=n_undistributed_samples, replace=True, p=perc_samples_per_cat)
    for cat_i in cats_to_distribute_samples :
        n_samples_per_cat[cat_i] += 1

    # Distribute samples of each cat to subjects.
    n_samples_per_cat_per_subj = np.zeros( [ n_sampl_categs, n_subjects_for_subep ] , dtype="int32" )
    for cat_i in range(n_sampl_categs) :
        n_samples_per_cat_per_subj[cat_i] += n_samples_per_cat[cat_i] // n_subjects_for_subep
        n_undistributed_samples_in_cat = n_samples_per_cat[cat_i] % n_subjects_for_subep
        for idx in range(n_undistributed_samples_in_cat):
            n_samples_per_cat_per_subj[cat_i, random.randint(0, n_subjects_for_subep - 1)] += 1

    return n_samples_per_cat_per_subj

def load_subj_and_get_samples(job_i,
                              log,
                              train_or_val,
                              run_input_checks,
                              cnn3d,
                              sampling_type,
                              paths_per_chan_per_subj,
                              paths_to_lbls_per_subj,
                              provided_mask,
                              paths_to_masks_per_subj,
                              provided_weightmaps_per_sampl_cat,
                              paths_to_wmaps_per_sampl_cat_per_subj,
                              # Pre-processing:
                              pad_input_imgs,
                              augm_params,

                              n_subjects_for_subep,
                              inds_of_subjects_for_subep,
                              n_samples_per_cat_per_subj
                              ):
    # paths_per_chan_per_subj: [ [ for channel-0 [ one path per subj ] ], ..., [ for channel-n  [ one path per subj ] ] ]
    # n_samples_per_cat_per_subj: np arr, shape [num sampling categories, num subjects in subepoch]
    # returns: ( channelsOfSegmentsFromThisSubjPerPathway, lblsOfPredictedPartOfSegmentsFromThisSubj )
    id_str = "[JOB:"+str(job_i)+"|PID:"+str(os.getpid())+"]"
    log.print3(id_str+" Load & sample from subject of index (in user's list): " + str(inds_of_subjects_for_subep[job_i]) +\
               " (Job #" + str(job_i) + "/" +str(n_subjects_for_subep)+")")

    # List, with [numberOfPathwaysThatTakeInput] sublists. Each sublist is list of [partImagesLoadedPerSubepoch] arrays [channels, R,C,Z].
    channelsOfSegmentsFromThisSubjPerPathway = [ [] for i in range(cnn3d.getNumPathwaysThatRequireInput()) ]
    lblsOfPredictedPartOfSegmentsFromThisSubj = [] # Labels only for the central/predicted part of segments.

    dimsOfPrimeSegmentRcz = cnn3d.pathways[0].getShapeOfInput(train_or_val)[2:]

    [allChannelsOfPatientInNpArray, # nparray(channels,dim0,dim1,dim2)
    gtLabelsImage,
    roiMask,
    arrayWithWeightMapsWhereToSampleForEachCategory,
    tupleOfPaddingPerAxesLeftRight #( (padLeftR, padRightR), (padLeftC,padRightC), (padLeftZ,padRightZ)). All 0s when no padding.
    ] = load_imgs_of_subject(log,
                             job_i,
                             train_or_val,
                             run_input_checks,
                             inds_of_subjects_for_subep[job_i],
                             paths_per_chan_per_subj,
                             True, # providedGtLabelsBool. If this getTheArr function is called (training), gtLabels should already been provided.
                             paths_to_lbls_per_subj,
                             cnn3d.num_classes,
                             provided_weightmaps_per_sampl_cat, # If true, must provide all. Placeholder in testing.
                             paths_to_wmaps_per_sampl_cat_per_subj, # Placeholder in testing.
                             provided_mask,
                             paths_to_masks_per_subj,
                             # Preprocessing
                             pad_input_imgs,
                             cnn3d.recFieldCnn, # used if pad_input_imgs
                             dimsOfPrimeSegmentRcz) # used if pad_input_imgs

    dimensionsOfImageChannel = allChannelsOfPatientInNpArray[0].shape
    finalWeightMapsToSampleFromPerCategoryForSubject = sampling_type.logicDecidingAndGivingFinalSamplingMapsForEachCategory(
                                                                                provided_weightmaps_per_sampl_cat,
                                                                                arrayWithWeightMapsWhereToSampleForEachCategory,
                                                                                True, #providedGtLabelsBool. True both for train and val.
                                                                                gtLabelsImage,
                                                                                provided_mask,
                                                                                roiMask,
                                                                                dimensionsOfImageChannel)
    str_samples_per_cat = " Got samples per category: "
    for cat_i in range( sampling_type.getNumberOfCategoriesToSample() ) :
        catString = sampling_type.getStringsPerCategoryToSample()[cat_i]
        numOfSegmsToExtractForThisCatFromThisSubject = n_samples_per_cat_per_subj[cat_i][job_i]
        finalWeightMapToSampleFromForThisCat = finalWeightMapsToSampleFromPerCategoryForSubject[cat_i]

        # Check if the weight map is fully-zeros. In this case, don't call the sampling function, just continue.
        # Note that this way, the data loaded on GPU will not be as much as I initially wanted. Thus calculate number-of-batches from this actual number of extracted segments.
        if np.sum(finalWeightMapToSampleFromForThisCat>0) == 0 :
            log.print3(id_str+" WARN: Sampling mask/map was found just zeros! No [" + catString + "] samples for this subject!")
            continue

        coords_of_samples = sample_coords_of_segments(log,
                                            job_i,
                                            numOfSegmsToExtractForThisCatFromThisSubject,
                                            dimsOfPrimeSegmentRcz,
                                            dimensionsOfImageChannel,
                                            finalWeightMapToSampleFromForThisCat)
        str_samples_per_cat += "[" + catString + ": " + str(len(coords_of_samples[0][0])) + "/" + str(numOfSegmsToExtractForThisCatFromThisSubject) + "] "

        # Use the just sampled coordinates of slices to actually extract the segments (data) from the subject's images.
        for image_part_i in range(len(coords_of_samples[0][0])) :
            coordsOfCentralVoxelOfThisImPart = coords_of_samples[0][:,image_part_i]

            [channelsForSegmentPerPathway,
            lblsOfPredictedPartOfSegment # used to be gtLabelsForThisImagePart, before extracting only for the central voxels.
            ] = extractSegmentGivenSliceCoords(train_or_val,
                                               cnn3d,
                                               coordsOfCentralVoxelOfThisImPart,
                                               allChannelsOfPatientInNpArray,
                                               gtLabelsImage)

            # Augmentation of segments
            channelsForSegmentPerPathway, lblsOfPredictedPartOfSegment = augment_patch(channelsForSegmentPerPathway, lblsOfPredictedPartOfSegment, augm_params)

            for pathway_i in range(cnn3d.getNumPathwaysThatRequireInput()) :
                channelsOfSegmentsFromThisSubjPerPathway[pathway_i].append(channelsForSegmentPerPathway[pathway_i])
            lblsOfPredictedPartOfSegmentsFromThisSubj.append(lblsOfPredictedPartOfSegment)

    log.print3(id_str + str_samples_per_cat)
    return (channelsOfSegmentsFromThisSubjPerPathway, lblsOfPredictedPartOfSegmentsFromThisSubj)

# roi_mask_filename and roiMinusLesion_mask_filename can be passed "no". In this case, the corresponding return result is nothing.
# This is so because: the do_training() function only needs the roiMinusLesion_mask, whereas the do_testing() only needs the roi_mask.
def load_imgs_of_subject(log,
                         job_i, # None in testing.
                         train_val_or_test,
                         run_input_checks,
                         subj_i,
                         paths_per_chan_per_subj,
                         providedGtLabelsBool,
                         listOfFilepathsToGtLabelsOfEachPatient,
                         num_classes,
                         provided_weightmaps_per_sampl_cat, # If true, must provide all wmaps. Placeholder in testing.
                         paths_to_wmaps_per_sampl_cat_per_subj, # Placeholder in testing.
                         provided_mask,
                         paths_to_masks_per_subj,
                         # Preprocessing
                         pad_input_imgs,
                         cnnReceptiveField, # only used if pad_input_imgs
                         dimsOfPrimeSegmentRcz
                         ):
    #listOfNiiFilepathNames: should be a list of lists. Each sublist corresponds to one certain subject.
    #...Each sublist should have as many elements(strings-filenamePaths) as numberOfChannels, point to the channels of this patient.
    id_str = "[JOB:"+str(job_i)+"|PID:"+str(os.getpid())+"] " if job_i is not None else "" # is None in testing.

    if subj_i >= len(paths_per_chan_per_subj) :
        raise ValueError(id_str+" The argument 'subj_i' given is greater than the filenames given for the .nii folders!")

    log.print3(id_str+"Loading subject with 1st channel at: "+str(paths_per_chan_per_subj[subj_i][0]))

    numberOfNormalScaleChannels = len(paths_per_chan_per_subj[0])

    tupleOfPaddingPerAxesLeftRight = ((0,0), (0,0), (0,0)) #This will be given a proper value if padding is performed.

    if provided_mask :
        fullFilenamePathOfRoiMask = paths_to_masks_per_subj[subj_i]
        roiMask = loadVolume(fullFilenamePathOfRoiMask)

        [roiMask, tupleOfPaddingPerAxesLeftRight] = padCnnInputs(roiMask, cnnReceptiveField, dimsOfPrimeSegmentRcz) if pad_input_imgs else [roiMask, tupleOfPaddingPerAxesLeftRight]
    else :
        roiMask = None

    #Load the channels of the patient.
    niiDimensions = None
    allChannelsOfPatientInNpArray = None
    #The below has dimensions (channels, 2). Holds per channel: [value to add per voxel for mean norm, value to multiply for std renorm]
    howMuchToAddAndMultiplyForNormalizationAugmentationForEachChannel = np.ones( (numberOfNormalScaleChannels, 2), dtype="float32")
    for channel_i in range(numberOfNormalScaleChannels):
        fullFilenamePathOfChannel = paths_per_chan_per_subj[subj_i][channel_i]
        if fullFilenamePathOfChannel != "-" : #normal case, filepath was given.
            channelData = loadVolume(fullFilenamePathOfChannel)

            [channelData, tupleOfPaddingPerAxesLeftRight] = padCnnInputs(channelData, cnnReceptiveField, dimsOfPrimeSegmentRcz) if pad_input_imgs else [channelData, tupleOfPaddingPerAxesLeftRight]

            if not isinstance(allChannelsOfPatientInNpArray, (np.ndarray)) :
                #Initialize the array in which all the channels for the patient will be placed.
                niiDimensions = list(channelData.shape)
                allChannelsOfPatientInNpArray = np.zeros( (numberOfNormalScaleChannels, niiDimensions[0], niiDimensions[1], niiDimensions[2]))

            allChannelsOfPatientInNpArray[channel_i] = channelData
        else : # "-" was given in the config-listing file. Do Min-fill!
            log.print3(id_str+" WARN: No modality #"+str(channel_i)+" given. Will make input channel full of zeros.")
            allChannelsOfPatientInNpArray[channel_i] = -4.0

    #Load the class labels.
    if providedGtLabelsBool : #For training (exact target labels) or validation on samples labels.
        fullFilenamePathOfGtLabels = listOfFilepathsToGtLabelsOfEachPatient[subj_i]
        imageGtLabels = loadVolume(fullFilenamePathOfGtLabels)

        if imageGtLabels.dtype.kind not in ['i','u']:
            log.print3(id_str+" WARN: Loaded labels are dtype ["+str(imageGtLabels.dtype)+"]. Rounding and casting them to int!")
            imageGtLabels = np.rint(imageGtLabels).astype("int32")

        if run_input_checks:
            check_gt_vs_num_classes(log, imageGtLabels, num_classes)

        [imageGtLabels, tupleOfPaddingPerAxesLeftRight] = padCnnInputs(imageGtLabels, cnnReceptiveField, dimsOfPrimeSegmentRcz) if pad_input_imgs else [imageGtLabels, tupleOfPaddingPerAxesLeftRight]
    else :
        imageGtLabels = None #For validation and testing

    if train_val_or_test != "test" and provided_weightmaps_per_sampl_cat==True : # in testing these weightedMaps are never provided, they are for training/validation only.
        n_sampl_categs = len(paths_to_wmaps_per_sampl_cat_per_subj)
        arrayWithWeightMapsWhereToSampleForEachCategory = np.zeros( [n_sampl_categs] + list(allChannelsOfPatientInNpArray[0].shape), dtype="float32" )
        for cat_i in range( n_sampl_categs ) :
            filepathsToTheWeightMapsOfAllPatientsForThisCategory = paths_to_wmaps_per_sampl_cat_per_subj[cat_i]
            filepathToTheWeightMapOfThisPatientForThisCategory = filepathsToTheWeightMapsOfAllPatientsForThisCategory[subj_i]
            weightedMapForThisCatData = loadVolume(filepathToTheWeightMapOfThisPatientForThisCategory)

            [weightedMapForThisCatData, tupleOfPaddingPerAxesLeftRight] = padCnnInputs(weightedMapForThisCatData, cnnReceptiveField, dimsOfPrimeSegmentRcz) if pad_input_imgs else [weightedMapForThisCatData, tupleOfPaddingPerAxesLeftRight]

            arrayWithWeightMapsWhereToSampleForEachCategory[cat_i] = weightedMapForThisCatData
    else :
        arrayWithWeightMapsWhereToSampleForEachCategory = None

    return [allChannelsOfPatientInNpArray, imageGtLabels, roiMask, arrayWithWeightMapsWhereToSampleForEachCategory, tupleOfPaddingPerAxesLeftRight]

#made for 3d
def sample_coords_of_segments(  log,
                                job_i,
                                numOfSegmentsToExtractForThisSubject,
                                dimsOfSegmentRcz,
                                dimensionsOfImageChannel,
                                weightMapToSampleFrom ) :
    """
    This function returns the coordinates (index) of the "central" voxel of sampled image parts (1voxel to the left if even part-dimension).
    It also returns the indices of the image parts, left and right indices, INCLUSIVE BOTH SIDES.

    Return value: [ rcz-coordsOfCentralVoxelsOfPartsSampled, rcz-sliceCoordsOfImagePartsSampled ]
    > coordsOfCentralVoxelsOfPartsSampled : an array with shape: 3(rcz) x numOfSegmentsToExtractForThisSubject.
        Example: [ rCoordsForCentralVoxelOfEachPart, cCoordsForCentralVoxelOfEachPart, zCoordsForCentralVoxelOfEachPart ]
        >> r/c/z-CoordsForCentralVoxelOfEachPart : A 1-dim array with numOfSegmentsToExtractForThisSubject, that holds the r-index within the image of each sampled part.
    > sliceCoordsOfImagePartsSampled : 3(rcz) x NumberOfImagePartSamples x 2. The last dimension has [0] for the lower boundary of the slice, and [1] for the higher boundary. INCLUSIVE BOTH SIDES.
        Example: [ r-sliceCoordsOfImagePart, c-sliceCoordsOfImagePart, z-sliceCoordsOfImagePart ]
    """
    id_str = "[JOB:"+str(job_i)+"|PID:"+str(os.getpid())+"]"
    # Check if the weight map is fully-zeros. In this case, return no element.
    # Note: Currently, the caller function is checking this case already and does not let this being called. Which is still fine.
    if np.sum(weightMapToSampleFrom>0) == 0 :
        log.print3(id_str+" WARN: The sampling mask/map was found just zeros! No image parts were sampled for this subject!")
        return [ [[],[],[]], [[],[],[]] ]

    imagePartsSampled = []

    #Now out of these, I need to randomly select one, which will be an ImagePart's central voxel.
    #But I need to be CAREFUL and get one that IS NOT closer to the image boundaries than the dimensions of the ImagePart permit.

    #I look for lesions that are not closer to the image boundaries than the ImagePart dimensions allow.
    #KernelDim is always odd. BUT ImagePart dimensions can be odd or even.
    #If odd, ok, floor(dim/2) from central.
    #If even, dim/2-1 voxels towards the begining of the axis and dim/2 towards the end. Ie, "central" imagePart voxel is 1 closer to begining.
    #BTW imagePartDim takes kernel into account (ie if I want 9^3 voxels classified per imagePart with kernel 5x5, I want 13 dim ImagePart)

    halfImagePartBoundaries = np.zeros( (len(dimsOfSegmentRcz), 2) , dtype='int32') #dim1: 1 row per r,c,z. Dim2: left/right width not to sample from (=half segment).

    #The below starts all zero. Will be Multiplied by other true-false arrays expressing if the relevant voxels are within boundaries.
    #In the end, the final vector will be true only for the indices of lesions that are within all boundaries.
    booleanNpArray_voxelsToCentraliseImPartsWithinBoundaries = np.zeros(weightMapToSampleFrom.shape, dtype="int32")

    # This loop leads to booleanNpArray_voxelsToCentraliseImPartsWithinBoundaries to be true for the indices ...
    # ...that allow getting an imagePart CENTERED on them, and be safely within image boundaries. Note that ...
    # ... if the imagePart is of even dimension, the "central" voxel is one voxel to the left.
    for rcz_i in range( len(dimsOfSegmentRcz) ) :
        if dimsOfSegmentRcz[rcz_i]%2 == 0: #even
            dimensionDividedByTwo = dimsOfSegmentRcz[rcz_i]//2
            halfImagePartBoundaries[rcz_i] = [dimensionDividedByTwo - 1, dimensionDividedByTwo] #central of ImagePart is 1 vox closer to begining of axes.
        else: #odd
            dimensionDividedByTwoFloor = math.floor(dimsOfSegmentRcz[rcz_i]//2) #eg 5/2 = 2, with the 3rd voxel being the "central"
            halfImagePartBoundaries[rcz_i] = [dimensionDividedByTwoFloor, dimensionDividedByTwoFloor]
    #used to be [halfImagePartBoundaries[0][0]: -halfImagePartBoundaries[0][1]], but in 2D case halfImagePartBoundaries might be ==0, causes problem and you get a null slice.
    booleanNpArray_voxelsToCentraliseImPartsWithinBoundaries[halfImagePartBoundaries[0][0]: dimensionsOfImageChannel[0] - halfImagePartBoundaries[0][1],
                                                            halfImagePartBoundaries[1][0]: dimensionsOfImageChannel[1] - halfImagePartBoundaries[1][1],
                                                            halfImagePartBoundaries[2][0]: dimensionsOfImageChannel[2] - halfImagePartBoundaries[2][1]] = 1

    constrainedWithImageBoundariesMaskToSample = weightMapToSampleFrom * booleanNpArray_voxelsToCentraliseImPartsWithinBoundaries
    #normalize the probabilities to sum to 1, cause the function needs it as so.
    constrainedWithImageBoundariesMaskToSample = constrainedWithImageBoundariesMaskToSample / (1.0* np.sum(constrainedWithImageBoundariesMaskToSample))

    flattenedConstrainedWithImageBoundariesMaskToSample = constrainedWithImageBoundariesMaskToSample.flatten()

    #This is going to be a 3xNumberOfImagePartSamples array.
    indicesInTheFlattenArrayThatWereSampledAsCentralVoxelsOfImageParts = np.random.choice(  constrainedWithImageBoundariesMaskToSample.size,
                                                                                size = numOfSegmentsToExtractForThisSubject,
                                                                                replace=True,
                                                                                p=flattenedConstrainedWithImageBoundariesMaskToSample)
    # np.unravel_index([listOfIndicesInFlattened], dims) returns a tuple of arrays (eg 3 of them if 3 dimImage),
    # where each of the array in the tuple has the same shape as the listOfIndices.
    # They have the r/c/z coords that correspond to the index of the flattened version.
    # So, coordsOfCentralVoxelsOfPartsSampled will be array of shape: 3(rcz) x numOfSegmentsToExtractForThisSubject.
    coordsOfCentralVoxelsOfPartsSampled = np.asarray(np.unravel_index(indicesInTheFlattenArrayThatWereSampledAsCentralVoxelsOfImageParts,
                                                                    constrainedWithImageBoundariesMaskToSample.shape #the shape of the roiMask/scan.
                                                                    )
                                                    )
    # Array with shape: 3(rcz) x NumberOfImagePartSamples x 2.
    # Last dimension has [0] for lowest boundary of slice, and [1] for highest boundary. INCLUSIVE BOTH SIDES.
    sliceCoordsOfImagePartsSampled = np.zeros(list(coordsOfCentralVoxelsOfPartsSampled.shape) + [2], dtype="int32")
    sliceCoordsOfImagePartsSampled[:,:,0] = coordsOfCentralVoxelsOfPartsSampled - halfImagePartBoundaries[ :, np.newaxis, 0 ] #np.newaxis broadcasts. To broadcast the -+.
    sliceCoordsOfImagePartsSampled[:,:,1] = coordsOfCentralVoxelsOfPartsSampled + halfImagePartBoundaries[ :, np.newaxis, 1 ]

    # coordsOfCentralVoxelsOfPartsSampled: Array of dimensions 3(rcz) x NumberOfImagePartSamples.
    # sliceCoordsOfImagePartsSampled: Array of dimensions 3(rcz) x NumberOfImagePartSamples x 2. ...
    # ... The last dim has [0] for the lower boundary of the slice, and [1] for the higher boundary.
    # ... The slice coordinates returned are INCLUSIVE BOTH sides.
    imagePartsSampled = [coordsOfCentralVoxelsOfPartsSampled, sliceCoordsOfImagePartsSampled]
    return imagePartsSampled

def getImagePartFromSubsampledImageForTraining( dimsOfPrimarySegment,
                                                recFieldCnn,
                                                subsampledImageChannels,
                                                image_part_slices_coords,
                                                subSamplingFactor,
                                                subsampledImagePartDimensions
                                                ) :
    """
    This returns an image part from the sampled data, given the image_part_slices_coords,
    which has the coordinates where the normal-scale image part starts and ends (inclusive).
    (Actually, in this case, the right (end) part of image_part_slices_coords is not used.)

    The way it works is NOT optimal. From the beginning of the normal-resolution part,
    it goes further to the left 1 receptive-field and then forward xSubsamplingFactor receptive-fields.
    This stops it from being used with arbitrary size of subsampled segment (decoupled by the high-res segment).
    Now, the subsampled patch has to be of the same size as the normal-scale.
    To change this, I should find where THE FIRST TOP LEFT CENTRAL (predicted) VOXEL is,
    and do the back-one-(sub)patch + front-3-(sub)patches from there, not from the begining of the patch.

    Current way it works (correct):
    If I have eg subsample factor=3 and 9 central-pred-voxels, I get 3 "central" voxels/patches for the
    subsampled-part. Straightforward. If I have a number of central voxels that is not an exact multiple of
    the subfactor, eg 10 central-voxels, I get 3+1 central voxels in the subsampled-part.
    When the cnn is convolving them, they will get repeated to 4(last-layer-neurons)*3(factor) = 12,
    and will get sliced down to 10, in order to have same dimension with the 1st pathway.
    """
    subsampledImageDimensions = subsampledImageChannels[0].shape

    subsampledChannelsForThisImagePart = np.ones((len(subsampledImageChannels),
                                                  subsampledImagePartDimensions[0],
                                                  subsampledImagePartDimensions[1],
                                                  subsampledImagePartDimensions[2]),
                                                 dtype = 'float32')

    numberOfCentralVoxelsClassifiedForEachImagePart_rDim = dimsOfPrimarySegment[0] - recFieldCnn[0] + 1
    numberOfCentralVoxelsClassifiedForEachImagePart_cDim = dimsOfPrimarySegment[1] - recFieldCnn[1] + 1
    numberOfCentralVoxelsClassifiedForEachImagePart_zDim = dimsOfPrimarySegment[2] - recFieldCnn[2] + 1

    #Calculate the slice that I should get, and where I should put it in the imagePart (eg if near the borders, and I cant grab a whole slice-imagePart).
    rSlotsPreviously = ((subSamplingFactor[0]-1)//2)*recFieldCnn[0] if subSamplingFactor[0]%2==1 \
                                                else (subSamplingFactor[0]-2)//2*recFieldCnn[0] + recFieldCnn[0]//2
    cSlotsPreviously = ((subSamplingFactor[1]-1)//2)*recFieldCnn[1] if subSamplingFactor[1]%2==1 \
                                                else (subSamplingFactor[1]-2)//2*recFieldCnn[1] + recFieldCnn[1]//2
    zSlotsPreviously = ((subSamplingFactor[2]-1)//2)*recFieldCnn[2] if subSamplingFactor[2]%2==1 \
                                                else (subSamplingFactor[2]-2)//2*recFieldCnn[2] + recFieldCnn[2]//2
    #1*17
    rToCentralVoxelOfAnAveragedArea = subSamplingFactor[0]//2 if subSamplingFactor[0]%2==1 else (subSamplingFactor[0]//2 - 1) #one closer to the beginning of dim. Same happens when I get parts of image.
    cToCentralVoxelOfAnAveragedArea = subSamplingFactor[1]//2 if subSamplingFactor[1]%2==1 else (subSamplingFactor[1]//2 - 1)
    zToCentralVoxelOfAnAveragedArea =  subSamplingFactor[2]//2 if subSamplingFactor[2]%2==1 else (subSamplingFactor[2]//2 - 1)
    #This is where to start taking voxels from the subsampled image. From the beginning of the imagePart(1 st patch)...
    #... go forward a few steps to the voxel that is like the "central" in this subsampled (eg 3x3) area.
    #...Then go backwards -Patchsize to find the first voxel of the subsampled.
    rlow = image_part_slices_coords[0][0] + rToCentralVoxelOfAnAveragedArea - rSlotsPreviously#These indices can run out of image boundaries. I ll correct them afterwards.
    #If the patch is 17x17, I want a 17x17 subsampled Patch. BUT if the imgPART is 25x25 (9voxClass), I want 3 subsampledPatches in my subsampPart to cover this area!
    #That is what the last term below is taking care of.
    #CAST TO INT because ceil returns a float, and later on when computing rHighNonInclToPutTheNotPaddedInSubsampledImPart I need to do INTEGER DIVISION.
    rhighNonIncl = int(rlow + subSamplingFactor[0]*recFieldCnn[0] + (math.ceil((numberOfCentralVoxelsClassifiedForEachImagePart_rDim*1.0)/subSamplingFactor[0]) - 1) * subSamplingFactor[0]) # excluding index in segment
    clow = image_part_slices_coords[1][0] + cToCentralVoxelOfAnAveragedArea - cSlotsPreviously
    chighNonIncl = int(clow + subSamplingFactor[1]*recFieldCnn[1] + (math.ceil((numberOfCentralVoxelsClassifiedForEachImagePart_cDim*1.0)/subSamplingFactor[1]) - 1) * subSamplingFactor[1])
    zlow = image_part_slices_coords[2][0] + zToCentralVoxelOfAnAveragedArea - zSlotsPreviously
    zhighNonIncl = int(zlow + subSamplingFactor[2]*recFieldCnn[2] + (math.ceil((numberOfCentralVoxelsClassifiedForEachImagePart_zDim*1.0)/subSamplingFactor[2]) - 1) * subSamplingFactor[2])

    rlowCorrected = max(rlow, 0)
    clowCorrected = max(clow, 0)
    zlowCorrected = max(zlow, 0)
    rhighNonInclCorrected = min(rhighNonIncl, subsampledImageDimensions[0])
    chighNonInclCorrected = min(chighNonIncl, subsampledImageDimensions[1])
    zhighNonInclCorrected = min(zhighNonIncl, subsampledImageDimensions[2]) #This gave 7

    rLowToPutTheNotPaddedInSubsampledImPart = 0 if rlow >= 0 else abs(rlow)//subSamplingFactor[0]
    cLowToPutTheNotPaddedInSubsampledImPart = 0 if clow >= 0 else abs(clow)//subSamplingFactor[1]
    zLowToPutTheNotPaddedInSubsampledImPart = 0 if zlow >= 0 else abs(zlow)//subSamplingFactor[2]

    dimensionsOfTheSliceOfSubsampledImageNotPadded = [  int(math.ceil((rhighNonInclCorrected - rlowCorrected)*1.0/subSamplingFactor[0])),
                                                        int(math.ceil((chighNonInclCorrected - clowCorrected)*1.0/subSamplingFactor[1])),
                                                        int(math.ceil((zhighNonInclCorrected - zlowCorrected)*1.0/subSamplingFactor[2]))
                                                        ]

    #I now have exactly where to get the slice from and where to put it in the new array.
    for channel_i in range(len(subsampledImageChannels)) :
        intensityZeroOfChannel = calculateTheZeroIntensityOf3dImage(subsampledImageChannels[channel_i])
        subsampledChannelsForThisImagePart[channel_i] *= intensityZeroOfChannel

        sliceOfSubsampledImageNotPadded = subsampledImageChannels[channel_i][   rlowCorrected : rhighNonInclCorrected : subSamplingFactor[0],
                                                                                clowCorrected : chighNonInclCorrected : subSamplingFactor[1],
                                                                                zlowCorrected : zhighNonInclCorrected : subSamplingFactor[2]
                                                                            ]
        subsampledChannelsForThisImagePart[
            channel_i,
            rLowToPutTheNotPaddedInSubsampledImPart : rLowToPutTheNotPaddedInSubsampledImPart+dimensionsOfTheSliceOfSubsampledImageNotPadded[0],
            cLowToPutTheNotPaddedInSubsampledImPart : cLowToPutTheNotPaddedInSubsampledImPart+dimensionsOfTheSliceOfSubsampledImageNotPadded[1],
            zLowToPutTheNotPaddedInSubsampledImPart : zLowToPutTheNotPaddedInSubsampledImPart+dimensionsOfTheSliceOfSubsampledImageNotPadded[2]] = sliceOfSubsampledImageNotPadded

    #placeholderReturn = np.ones([3,19,19,19], dtype="float32") #channel, dims
    return subsampledChannelsForThisImagePart

def shuffleSegmentsOfThisSubepoch(  channelsOfSegmentsForSubepochListPerPathway,
                                    lblsForPredictedPartOfSegmentsForSubepochList ) :
    numOfPathwayWithInput = len(channelsOfSegmentsForSubepochListPerPathway)
    inputToZip = [ sublistForPathway for sublistForPathway in channelsOfSegmentsForSubepochListPerPathway ]
    inputToZip += [ lblsForPredictedPartOfSegmentsForSubepochList ]

    combined = list(zip(*inputToZip)) #list() for python3 compatibility, as range cannot get assignment in shuffle()
    random.shuffle(combined)
    shuffledInputListsToZip = list(zip(*combined))

    shuffledChannelsOfSegmentsForSubepochListPerPathway = [ sublistForPathway for sublistForPathway in shuffledInputListsToZip[:numOfPathwayWithInput] ]
    shuffledLblsForPredictedPartOfSegmentsForSubepochList = shuffledInputListsToZip[numOfPathwayWithInput]

    return [shuffledChannelsOfSegmentsForSubepochListPerPathway, shuffledLblsForPredictedPartOfSegmentsForSubepochList]

# I must merge this with function: extractSegmentsGivenSliceCoords() that is used for Testing! Should be easy!
# This is used in training/val only.
def extractSegmentGivenSliceCoords(train_or_val,
                                   cnn3d,
                                   coordsOfCentralVoxelOfThisImPart,
                                   allChannelsOfPatientInNpArray,
                                   gtLabelsImage) :
    # allChannelsOfPatientInNpArray: numpy array [ n_channels, x, y, z ]

    channelsForSegmentPerPathway = []
    # Sampling
    for pathway in cnn3d.pathways[:1] : #Hack. The rest of this loop can work for the whole .pathways...
        # ... BUT the loop does not check what happens if boundaries are out of limits, to fill with zeros. This is done in getImagePartFromSubsampledImageForTraining().
        #... Update it in a nice way to be done here, and then take getImagePartFromSubsampledImageForTraining() out and make loop go for every pathway.

        if pathway.pType() == pt.FC :
            continue
        subSamplingFactor = pathway.subsFactor()
        pathwayInputShapeRcz = pathway.getShapeOfInput(train_or_val)[2:]
        leftBoundaryRcz = [ coordsOfCentralVoxelOfThisImPart[0] - subSamplingFactor[0]*(pathwayInputShapeRcz[0]-1)//2,
                            coordsOfCentralVoxelOfThisImPart[1] - subSamplingFactor[1]*(pathwayInputShapeRcz[1]-1)//2,
                            coordsOfCentralVoxelOfThisImPart[2] - subSamplingFactor[2]*(pathwayInputShapeRcz[2]-1)//2]
        rightBoundaryRcz = [leftBoundaryRcz[0] + subSamplingFactor[0]*pathwayInputShapeRcz[0],
                            leftBoundaryRcz[1] + subSamplingFactor[1]*pathwayInputShapeRcz[1],
                            leftBoundaryRcz[2] + subSamplingFactor[2]*pathwayInputShapeRcz[2]]

        channelsForThisImagePart = allChannelsOfPatientInNpArray[:,
                                                                leftBoundaryRcz[0] : rightBoundaryRcz[0] : subSamplingFactor[0],
                                                                leftBoundaryRcz[1] : rightBoundaryRcz[1] : subSamplingFactor[1],
                                                                leftBoundaryRcz[2] : rightBoundaryRcz[2] : subSamplingFactor[2]]

        channelsForSegmentPerPathway.append(channelsForThisImagePart)

    # Extract the samples for secondary pathways. This whole for can go away, if I update above code to check to slices out of limits.
    for pathway_i in range(len(cnn3d.pathways)) : # Except Normal 1st, cause that was done already.
        if cnn3d.pathways[pathway_i].pType() == pt.FC or cnn3d.pathways[pathway_i].pType() == pt.NORM:
            continue
        #this datastructure is similar to channelsForThisImagePart, but contains voxels from the subsampled image.
        dimsOfPrimarySegment = cnn3d.pathways[pathway_i].getShapeOfInput(train_or_val)[2:]
        slicesCoordsOfSegmForPrimaryPathway = [ [leftBoundaryRcz[0], rightBoundaryRcz[0]-1], [leftBoundaryRcz[1], rightBoundaryRcz[1]-1], [leftBoundaryRcz[2], rightBoundaryRcz[2]-1] ] # rightmost  are placeholders here.
        channsForThisSubsampledPartAndPathway = getImagePartFromSubsampledImageForTraining(dimsOfPrimarySegment=dimsOfPrimarySegment,
                                                                                recFieldCnn=cnn3d.recFieldCnn,
                                                                                subsampledImageChannels=allChannelsOfPatientInNpArray,
                                                                                image_part_slices_coords=slicesCoordsOfSegmForPrimaryPathway,
                                                                                subSamplingFactor=cnn3d.pathways[pathway_i].subsFactor(),
                                                                                subsampledImagePartDimensions=cnn3d.pathways[pathway_i].getShapeOfInput(train_or_val)[2:]
                                                                                )

        channelsForSegmentPerPathway.append(channsForThisSubsampledPartAndPathway)

    # Get ground truth labels for training.
    numOfCentralVoxelsClassifRcz = cnn3d.finalTargetLayer_outputShape[train_or_val][2:]
    leftBoundaryRcz = [ coordsOfCentralVoxelOfThisImPart[0] - (numOfCentralVoxelsClassifRcz[0]-1)//2,
                        coordsOfCentralVoxelOfThisImPart[1] - (numOfCentralVoxelsClassifRcz[1]-1)//2,
                        coordsOfCentralVoxelOfThisImPart[2] - (numOfCentralVoxelsClassifRcz[2]-1)//2]
    rightBoundaryRcz = [leftBoundaryRcz[0] + numOfCentralVoxelsClassifRcz[0],
                        leftBoundaryRcz[1] + numOfCentralVoxelsClassifRcz[1],
                        leftBoundaryRcz[2] + numOfCentralVoxelsClassifRcz[2]]
    lblsOfPredictedPartOfSegment = gtLabelsImage[ leftBoundaryRcz[0] : rightBoundaryRcz[0],
                                                  leftBoundaryRcz[1] : rightBoundaryRcz[1],
                                                  leftBoundaryRcz[2] : rightBoundaryRcz[2] ]

    # Make COPIES of the segments, instead of having a VIEW (slice) of them. This is so that the the whole volume are afterwards released from RAM.
    channelsForSegmentPerPathway = [ np.array(pathw_channs, copy=True, dtype='float32') for pathw_channs in channelsForSegmentPerPathway  ]
    lblsOfPredictedPartOfSegment = np.copy(lblsOfPredictedPartOfSegment)

    return [ channelsForSegmentPerPathway, lblsOfPredictedPartOfSegment ]

#################################################################################################################################
#                                                                                #
#       Below are functions for testing only. There is duplication with training. They are not the same, but could be merged.   #
#                                                                                #
#################################################################################################################################

# This is very similar to sample_coords_of_segments() I believe, which is used for training. Consider way to merge them.
def getCoordsOfAllSegmentsOfAnImage(log,
                                    dimsOfPrimarySegment, # RCZ dims of input to primary pathway (NORMAL). Which should be the first one in .pathways.
                                    strideOfSegmentsPerDimInVoxels,
                                    batch_size,
                                    channelsOfImageNpArray,
                                    roiMask
                                    ) :
    # channelsOfImageNpArray: np array [n_channels, x, y, z]
    log.print3("Starting to (tile) extract Segments from the images of the subject for Segmentation...")

    sliceCoordsOfSegmentsToReturn = []

    niiDimensions = list(channelsOfImageNpArray[0].shape) # Dims of the volumes

    zLowBoundaryNext=0; zAxisCentralPartPredicted = False;
    while not zAxisCentralPartPredicted :
        zFarBoundary = min(zLowBoundaryNext+dimsOfPrimarySegment[2], niiDimensions[2]) #Excluding.
        zLowBoundary = zFarBoundary - dimsOfPrimarySegment[2]
        zLowBoundaryNext = zLowBoundaryNext + strideOfSegmentsPerDimInVoxels[2]
        zAxisCentralPartPredicted = False if zFarBoundary < niiDimensions[2] else True #THIS IS THE IMPORTANT CRITERION.

        cLowBoundaryNext=0; cAxisCentralPartPredicted = False;
        while not cAxisCentralPartPredicted :
            cFarBoundary = min(cLowBoundaryNext+dimsOfPrimarySegment[1], niiDimensions[1]) #Excluding.
            cLowBoundary = cFarBoundary - dimsOfPrimarySegment[1]
            cLowBoundaryNext = cLowBoundaryNext + strideOfSegmentsPerDimInVoxels[1]
            cAxisCentralPartPredicted = False if cFarBoundary < niiDimensions[1] else True

            rLowBoundaryNext=0; rAxisCentralPartPredicted = False;
            while not rAxisCentralPartPredicted :
                rFarBoundary = min(rLowBoundaryNext+dimsOfPrimarySegment[0], niiDimensions[0]) #Excluding.
                rLowBoundary = rFarBoundary - dimsOfPrimarySegment[0]
                rLowBoundaryNext = rLowBoundaryNext + strideOfSegmentsPerDimInVoxels[0]
                rAxisCentralPartPredicted = False if rFarBoundary < niiDimensions[0] else True

                if isinstance(roiMask, (np.ndarray)) : #In case I pass a brain-mask, I ll use it to only predict inside it. Otherwise, whole image.
                    if not np.any(roiMask[rLowBoundary:rFarBoundary,
                                            cLowBoundary:cFarBoundary,
                                            zLowBoundary:zFarBoundary
                                            ]) : #all of it is out of the brain so skip it.
                        continue

                sliceCoordsOfSegmentsToReturn.append([ [rLowBoundary, rFarBoundary-1], [cLowBoundary, cFarBoundary-1], [zLowBoundary, zFarBoundary-1] ])

    #I need to have a total number of image-parts that can be exactly-divided by the 'batch_size'. For this reason, I add in the far end of the list multiple copies of the last element.
    total_number_of_image_parts = len(sliceCoordsOfSegmentsToReturn)
    number_of_imageParts_missing_for_exact_division =  batch_size - total_number_of_image_parts%batch_size if total_number_of_image_parts%batch_size != 0 else 0
    for extra_useless_image_part_i in range(number_of_imageParts_missing_for_exact_division) :
        sliceCoordsOfSegmentsToReturn.append(sliceCoordsOfSegmentsToReturn[-1])

    #I think that since the parts are acquired in a certain order and are sorted this way in the list, it is easy
    #to know which part of the image they came from, as it depends only on the stride-size and the imagePart size.

    log.print3("Finished (tiling) extracting Segments from the images of the subject for Segmentation.")

    # sliceCoordsOfSegmentsToReturn: list with 3 dimensions. numberOfSegments x 3(rcz) x 2 (lower and upper limit of the segment, INCLUSIVE both sides)
    return [sliceCoordsOfSegmentsToReturn]

# I must merge this with function: extractSegmentGivenSliceCoords() that is used for Training/Validation! Should be easy!
# This is used in testing only.
def extractSegmentsGivenSliceCoords(cnn3d,
                                    sliceCoordsOfSegmentsToExtract,
                                    channelsOfImageNpArray,
                                    recFieldCnn) :
    # channelsOfImageNpArray: numpy array [ n_channels, x, y, z ]
    numberOfSegmentsToExtract = len(sliceCoordsOfSegmentsToExtract)
    channsForSegmentsPerPathToReturn = [ [] for i in range(cnn3d.getNumPathwaysThatRequireInput()) ] # [pathway, image parts, channels, r, c, z]
    dimsOfPrimarySegment = cnn3d.pathways[0].getShapeOfInput("test")[2:] # RCZ dims of input to primary pathway (NORMAL). Which should be the first one in .pathways.

    for segment_i in range(numberOfSegmentsToExtract) :
        rLowBoundary = sliceCoordsOfSegmentsToExtract[segment_i][0][0]; rFarBoundary = sliceCoordsOfSegmentsToExtract[segment_i][0][1]
        cLowBoundary = sliceCoordsOfSegmentsToExtract[segment_i][1][0]; cFarBoundary = sliceCoordsOfSegmentsToExtract[segment_i][1][1]
        zLowBoundary = sliceCoordsOfSegmentsToExtract[segment_i][2][0]; zFarBoundary = sliceCoordsOfSegmentsToExtract[segment_i][2][1]
        # segment for primary pathway
        channsForPrimaryPathForThisSegm = channelsOfImageNpArray[:,
                                                                rLowBoundary:rFarBoundary+1,
                                                                cLowBoundary:cFarBoundary+1,
                                                                zLowBoundary:zFarBoundary+1
                                                                ]
        channsForSegmentsPerPathToReturn[0].append(channsForPrimaryPathForThisSegm)

        #Subsampled pathways
        for pathway_i in range(len(cnn3d.pathways)) : # Except Normal 1st, cause that was done already.
            if cnn3d.pathways[pathway_i].pType() == pt.FC or cnn3d.pathways[pathway_i].pType() == pt.NORM:
                continue
            slicesCoordsOfSegmForPrimaryPathway = [ [rLowBoundary, rFarBoundary-1], [cLowBoundary, cFarBoundary-1], [zLowBoundary, zFarBoundary-1] ] #the right hand values are placeholders in this case.
            channsForThisSubsPathForThisSegm = getImagePartFromSubsampledImageForTraining(  dimsOfPrimarySegment=dimsOfPrimarySegment,
                                                                                recFieldCnn=recFieldCnn,
                                                                                subsampledImageChannels=channelsOfImageNpArray,
                                                                                image_part_slices_coords=slicesCoordsOfSegmForPrimaryPathway,
                                                                                subSamplingFactor=cnn3d.pathways[pathway_i].subsFactor(),
                                                                                subsampledImagePartDimensions=cnn3d.pathways[pathway_i].getShapeOfInput("test")[2:]
                                                                                )
            channsForSegmentsPerPathToReturn[pathway_i].append(channsForThisSubsPathForThisSegm)

    return [channsForSegmentsPerPathToReturn]

###########################################
# Checks whether the data is as expected  #
###########################################

def check_gt_vs_num_classes(log, img_gt, num_classes):
    id_str = "["+str(os.getpid())+"]"
    max_in_gt = np.max(img_gt)
    if np.max(img_gt) > num_classes-1 : # num_classes includes background=0
        msg =  id_str+" ERROR:\t GT labels included a label value ["+str(max_in_gt)+"] greater than what CNN expects."+\
                "\n\t In model-config the number of classes was specified as ["+str(num_classes)+"]."+\
                "\n\t Check your data or change the number of classes in model-config."+\
                "\n\t Note: number of classes in model config should include the background as a class."
        log.print3(msg)
        raise ValueError(msg)

max_stars_repo_name: YongLiuLab/BrainRadiomicsTools

max_issues_repo_name: YongLiuLab/BrainRadiomicsTools

max_forks_repo_name: YongLiuLab/BrainRadiomicsTools


text:
"""Generate data with bandit feedback."""
import copy
from typing import List, Optional, Tuple, Union

import numpy as np
from pandas import DataFrame

from pymab.bandit import BaseBandit
from pymab.policy import BaseContextualPolicy, BasePolicy, BaseThompsonSampling

class DataGenerator():
    """Data generator.

    Parameters
    ----------
    n_arms: int
        The number of given bandit arms.

    n_features: int
        The dimention of context vectors.

    noise: int, optional(default=None)
        The variance of gaussian noise on linear model of contextual rewards.

    n_rounds: int
        The number of rounds in a simulation.

    randomized: bool, optional(default=False)
        Whether the past policy is randomized or not.

    """

    def __init__(self,
                 policy: Union[BasePolicy, BaseThompsonSampling, BaseContextualPolicy],
                 bandit: BaseBandit,
                 n_rounds: int,
                 randomized: bool=False) -> None:
        """Initialize Class."""
        self.policy = policy
        self.bandit = bandit
        self.n_arms = self.bandit.n_arms
        self.n_features = self.bandit.n_features
        self.n_rounds = n_rounds
        self.randomized = randomized

    def generate_data(self) -> Tuple:
        """Generate data."""
        data = np.concatenate([np.expand_dims(np.random.randint(2, size=self.n_features), axis=1)
                               for i in np.arange(self.n_rounds)], axis=1)

        chosen_arms = np.zeros(self.n_rounds, dtype=int)
        rewards = np.zeros(self.n_rounds)
        p, b = copy.deepcopy(self.policy), copy.deepcopy(self.bandit)

        for t in np.arange(self.n_rounds):
            x = np.copy(data[:, t])

            if self.randomized:
                chosen_arm = np.random.randint(self.n_arms)
            else:
                chosen_arm = p.select_arm(x) if p._policy_type == "contextual" else p.select_arm()

            reward = b.pull(x=x, chosen_arm=chosen_arm)
            chosen_arms[t], rewards[t] = chosen_arm, reward
            p.update(x, chosen_arm, reward) if p._policy_type == "contextual" else p.update(chosen_arm, reward)

        return data.T, chosen_arms, rewards

ext: py

max_stars_repo_name: smn-ailab/ysaito-qiita

max_issues_repo_name: smn-ailab/ysaito-qiita

max_forks_repo_name: smn-ailab/ysaito-qiita

reason: import numpy


text:
# Processing (proc) module for neural analysis package

# Load MGR data
def load_mgr_data():
    # Import libraries
    import pandas as pd

    # Define file path and load
    file_path = '/Volumes/Samsung_T5/Random Datasets/Ike_PMd_MemoryGuidedReach/Ike_MGR.hdf'
    df = pd.read_hdf(file_path)

    return df

# Calculate movement onset
def calc_movement_onset(df):
    import numpy as np
    import plottools as pt
    import matplotlib as mpl
    mpl.use('PDF')
    import matplotlib.pyplot as plt
    import copy

    # Define speed threshold used for movement onset
    s_thresh = 0.2

    # Create figure
    ax_hndl = pt.create_subplot(1, 1)

    # Iterate over trials and calculate speed
    n_trials = df.shape[0]
    s_all = []
    onset_idx = []
    for i in range(n_trials):
        # Get velocity, time, calculate speed
        v = df['vel'][i][0:2, :]
        t = df['time'][i]
        s = np.sqrt(np.sum(v**2, 0))

        # First get timing events -- trajectory onset/offset
        traj_onset = df['trajectoryOnset'][i]
        traj_offset = df['trajectoryOffset'][i]

        # Get time of max speed.  Find this between trajectory onset and offset
        s_temp = copy.copy(s)  # Create a copy of s
        s_temp[t < traj_onset] = 0  # zero out all speeds pre-trajectory onset (movement onset can't occur here)
        s_temp[t > traj_offset] = 0  # also zero-out all speeds post-trajectory offset
        max_ind = np.argmax(s_temp)  # Index corresponding to the max speed
        s_max = s[max_ind]
        t_max = t[max_ind]

        # Define masks used to find movement onset
        s_mask = s < (s_max * s_thresh)  # All speeds less than the threshold
        t_onset_mask = t > traj_onset
        t_max_mask = t < t_max
        mask = s_mask & t_onset_mask & t_max_mask  # All time points after trajectory onset but before max speed
        valid_idx = np.nonzero(mask)  # Get list of valid indices (speeds less than the threshold)
        onset_idx.append(valid_idx[0][-1] + 1)  # Movement onset will be the first time point after the last non-zero value

        # Add speed to list
        s_all.append(s)

    # Add movement onset and speed to dataframe
    df['speed'] = s_all
    df['onset_idx'] = onset_idx

    return None

# Bin neural and kinematic data
def bin_trial_data(
        df,
        bin_width = 20,  # Bin width (ms)
        bin_start = 0,  # Start of bin times (ms)
        valid_sort = range(1, 31),  # Valid ranges of sort codes to use
        kin_offset = 0  # Temporal offset to apply to kinematic data
        ):
    """
    Bin trial data.

    This function bins neural and kinematic data.

    """
    # To access a single element in the dataframe df, use 'df.iloc[idx]'
    # Iterate over trials in df, determine the list of valid channel/sort

    # Iterate over trials
    for trial in df.iterrows():
        # Iterate over channel/sort combinations.  Note that the 'iterrows()'
        # function is not the most efficient way to parse a dataframe.  The
        # 'apply()' method can be used to apply a specific function to each
        # element.
        for ch, srt, st in zip(
            trial.spike_channel,
            trial.spike_sort,
            trial.spike_times
        ):

            # Check to see if ch and srt are in the valid list.  If so, bin the
            # spike times.
            X = []

        # Bin hand position data as well
        Y = []

    # Outputs
    bin_data = []  # Binned neural spike counts (dict, one element per trial)

    return X, Y

# Bin data helper function
def bin_data(x, t, bin_width, bin_start):
    """
    Bin data.

    Bin the provided data into evenly-distributed time bins.

    """
    # Define time grid
    t_bins = []
    x_bins = []

    # Iterate through time bins

    return x_bins, t_bins

max_stars_repo_name: alandegenhart/decoding

max_issues_repo_name: alandegenhart/decoding

max_forks_repo_name: alandegenhart/decoding


text:
#!/usr/bin/env python
# import all needed libraries

import sys
import time
import sockjs.tornado
from tornado import web, ioloop
from sockjs.tornado import SockJSRouter, SockJSConnection
import pigpio
import subprocess
import os
import signal
from thread import start_new_thread
import random
import colorsys
import pyaudio
from scipy.signal import butter, lfilter, freqz
import numpy as np

# Initial setup of GPIO pins
pi = pigpio.pi()

############################### setting basic options  ###############################
bright    = 255

# The Pins. Use Broadcom numbers.
RED_PIN   = 17
GREEN_PIN = 22
BLUE_PIN  = 24

# listening port
port      = 1713

# Global variables for Music
CHUNK     = 512 # How many bytes of audio to read at a time
global HUE
HUE = 0

############################### setting effect options  ##############################
mode = "nothing"

############################### basic functions  ###############################
class LedController:
    def hex_to_rgb(self, hex):
        hex = hex.lstrip('#')
        lv = len(hex)
        rgb = tuple(int(hex[i:i + lv // 3], 16) for i in range(0, lv, lv // 3))
        return rgb

    def rgb_to_hex(self, rgb):
        hex = "#"
        for i in range(3):
            hex + hex(rgb[i]).split('x')[1]
        return hex

    def hsl_to_rgb(self, hsl):
        if any(i > 1 for i in hsl):
            hsl[0] /= 360
            hsl[1] /= 100
            hsl[2] /= 100
        colour_tuple = tuple(i * 255 for i in colorsys.hls_to_rgb(hsl[0], hsl[2], hsl[1]))
        return colour_tuple

    def checkRGB(self, color):
        if (color < 0):
            color = 0
        if (color > 255):
            color = 255
        return color

    def scale(self, brightness):
        realBrightness = int(int(brightness) * (float(bright) / 255.0))
        # Ensure we are giving correct values
        if realBrightness < 0:
            realBrightness = 0.0
        elif realBrightness > 100:
            realBrightness = 100.0
        return realBrightness

    def setLights(self, pin, brightness):
        realBrightness = self.scale(brightness)
        pi.set_PWM_dutycycle(pin, realBrightness)

    def setRGB(self, r, g, b):
        r = self.checkRGB(r)
        g = self.checkRGB(g)
        b = self.checkRGB(b)
        self.setLights(RED_PIN, r)
        self.setLights(GREEN_PIN, g)
        self.setLights(BLUE_PIN,b)
        print "changing color to rgb(" + str(r) + "," + str(g) + "," + str(b) + ")"
        #BrokerConnection.color_broadcaster(r,g,b)

    def noWhite(self, r, g, b):
        r /= 255.0
        g /= 255.0
        b /= 255.0
        HSL = colorsys.rgb_to_hls(r, g, b)
        h = HSL[0]
        l = HSL[1]
        s = HSL[2]
        if (l > 0.8):
            l *= 0.8 # scale down lightness when higher than 80%
        if (s < 0.4):
            s = (s * 0.6) + 0.4 # scale saturation up when lower than 40%
        return tuple(i * 255 for i in colorsys.hls_to_rgb(h, l, s))

############################### Rainbow functions  ###############################
class Rainbow:

    def updateColor(self, color, step):
            color += step
            if color > 255:
                color = 255
            if color < 0:
                color = 0
            return color

    def fader(self, r, g, b):
        lc = LedController()
        if not ((r == 255 or r == 0) and (b == 255 or b == 0) and (g == 255 or g == 0) and (r == 255 and g == 255 and b == 255) and (r == 0 and g == 0 and b == 0)):
            while r < 255:
                r = self.updateColor(r, STEPS)
                lc.setRGB(r, g, b)
            while b > 0:
                b = self.updateColor(b, -STEPS)
                lc.setRGB(r, g, b)
        while (mode == "Rainbow"):
            if r == 255 and b == 0 and g < 255:
                g = self.updateColor(g, STEPS)
                lc.setRGB(r, g, b)

            elif g == 255 and b == 0 and r > 0:
                r = self.updateColor(r, -STEPS)
                lc.setRGB(r, g, b)

            elif r == 0 and g == 255 and b < 255:
                b = self.updateColor(b, STEPS)
                lc.setRGB(r, g, b)

            elif r == 0 and b == 255 and g > 0:
                g = self.updateColor(g, -STEPS)
                lc.setRGB(r, g, b)

            elif g == 0 and b == 255 and r < 255:
                r = self.updateColor(r, STEPS)
                lc.setRGB(r, g, b)

            elif r == 255 and g == 0 and b > 0:
                b = self.updateColor(b, -STEPS)
                lc.setRGB(r, g, b)

############################### Music functions  ###############################

class FreqAnalyser:
    # Filtering based on
    # http://wiki.scipy.org/Cookbook/ButterworthBandpass

    def __init__(self, channels, sample_rate, leds=None):
        self.leds = leds # Not needed if just plotting
        self.channels = channels
        self.sample_rate = sample_rate
        self.nyquist = float(sample_rate) / 2

        # Filter order - higher the order the sharper
        # the curve
        order = 3

        # Cut off frequencies:
        # Low pass filter
        cutoff = 200 / self.nyquist
        # Numerator (b) and denominator (a)
        # polynomials of the filter.
        b, a = butter(order, cutoff, btype='lowpass')
        self.low_b = b
        self.low_a = a

        # High pass filter
        cutoff = 4000 / self.nyquist
        b, a = butter(order, cutoff, btype='highpass')
        self.high_b = b
        self.high_a = a

        # Keep track of max brightness for each
        # colour
        self.max = [0.0, 0.0, 0.0]
        # Make different frequencies fall faster
        # bass needs to be punchy.
        self.fall = [15.0, 2.5, 5.0]

    def filter(self, data):
        # Apply low filter
        self.low_data = lfilter(self.low_b,
                                self.low_a,
                                data)

        # Apply high filter
        self.high_data = lfilter(self.high_b,
                                 self.high_a,
                                 data)

        # Get mid data by doing signal - (low + high)
        self.mid_data = np.subtract(data,
                        np.add(self.low_data,
                               self.high_data))

    @staticmethod
    def rms(data):
        # Return root mean square of data set
        # (i.e. average amplitude)
        return np.sqrt(np.mean(np.square(data)))

    def change_leds(self):
        lc = LedController()

        # Get average amplitude
        l = []
        l.append(self.rms(self.low_data))
        l.append(self.rms(self.mid_data))
        l.append(self.rms(self.high_data))

        if mode == "Music":
            HUEcolor = MusicColor
            swift = (sum(l) * random.uniform(0, 7))
            if swift < 0.5:
                swift = 0.5

            if (HUEcolor == "Auto"):
                global HUE
                HUE += swift
            else:
                rgb = lc.hex_to_rgb(HUEcolor)
                rgb = [float(rgb[i]) / 255.0 for i in range(3)]
                global HUE
                HUE = colorsys.rgb_to_hls(rgb[0], rgb[1], rgb[2])[0] * 360
            light = 0.01 + l[0]
            if HUE > 360:
                HUE = 0 + (HUE - 360)
            if light > 0.6:
                light = 0.6

            RGB = lc.hsl_to_rgb([HUE / 360, 1, light])

        elif mode == "Music1":
            equalizer = MusicColor

            for i in range(3):
                # Do any number fudging to make it look better
                # here - probably want to avoid high values of
                # all because it will be white
                # (Emphasise/Reduce bass, mids, treble)
                l[i] *= float(equalizer[i])
                l[i] = (l[i] * 256) - 1

                # Use new val if > previous max
                if l[i] > self.max[i]:
                    self.max[i] = l[i]
                else:
                    # Otherwise, decrement max and use that
                    # Gives colour falling effect
                    self.max[i] -= self.fall[i]
                    if self.max[i] < 0:
                        self.max[i] = 0
                l[i] = self.max[i]
            RGB = l

        lc.setRGB(RGB[0], RGB[1], RGB[2])

class AudioController:
    def __init__(self, leds):
        self.line_in = True

        self.leds = leds
        self.p = pyaudio.PyAudio()

    def more(self):
        try:
            # Return line in data
            return self.stream.read(CHUNK)
        except:
            print "line-in error"
            return 'ab'

    def analyse(self, data):
        # Convert to numpy array and filter
        data = np.fromstring(data, dtype=np.int16)

        # Convert int16 to float for dsp
        data = np.float32(data/32768.0)

        # Send to filter
        self.analyser.filter(data)

        self.analyser.change_leds()

    def record_setup(self):
        self.channels = 1
        self.sample_rate = 44100
        self.stream = self.p.open(format = pyaudio.paInt16,
                                  channels = self.channels,
                                  rate = self.sample_rate,
                                  input=True,
                                  frames_per_buffer=CHUNK)

    def loop(self):
        # Main processing loop
        # Do appropriate setup
        self.record_setup()

        self.analyser = FreqAnalyser(self.channels,
                                     self.sample_rate,
                                     self.leds)

        # Read the first block of audio data
        data = self.more()

        # While there is still audio left
        while (mode == "Music") or (mode == "Music1"):
            try:
                # Analyse data and change LEDs
                self.analyse(data)

                # Get more audio data
                data = self.more()
            except KeyboardInterrupt:
                break

        # Tidy up
        self.stream.close()
        self.p.terminate()

############################### other Effects functions  ###############################

class Effects:
    def Flasher(self):
        while (mode == "Flash"):
            random = self.Random_color()
            RGB = LedController().noWhite(random[0],random[1],random[2])
            r = RGB[0]
            g = RGB[1]
            b = RGB[2]
            LedController().setRGB(r, g, b)
            time.sleep(tempo)

    def Strober(self):
        while (mode == "Strobe"):
            if tempo < 0:
                random = self.Random_color()
                LedController().setRGB(random[0],random[1],random[2])
                self.wait_s(tempo)
                LedController().setRGB(0, 0, 0)
                self.wait_s(tempo)
            else:
                LedController().setRGB(255, 255, 255)
                self.wait_s(tempo)
                LedController().setRGB(0, 0, 0)
                self.wait_s(tempo)

    def Random_color(self):
        h = random.uniform(0, 100) / 100
        s = random.uniform(95, 100) / 100
        v = random.uniform(88, 100) / 100
        return tuple(i * 255 for i in colorsys.hsv_to_rgb(h, s, v))

    def wait_s(self,seconds):
        if seconds < 0:
            time.sleep((seconds) * (-1))
        elif seconds >= 0:
            time.sleep(seconds)

###################################### Socket ######################################
class BrokerConnection(sockjs.tornado.SockJSConnection):
    clients = set()
    lc = LedController()
    rb = Rainbow()
    ef = Effects()
    ac = AudioController(lc)

    def on_open(self, info):
        # When new client comes in, will add it to the clients list
        self.clients.add(self)

    def on_message(self, message):
        # For every incoming message, broadcast it to all clients
        #self.broadcast(self.clients, message)
	    # Set RGB color
        if not 'r' in locals():
            r = 0
        if not 'g' in locals():
            g = 0
        if not 'b' in locals():
            b = 0
        aRGB = self.message_analyser(message)
        if len(aRGB) == 3:
            r = float(aRGB[0])
            g = float(aRGB[1])
            b = float(aRGB[2])
            if not mode == "nothing":
                mode = "nothing"
                time.sleep(0.1)
            self.lc.setRGB(r,g,b)

        if len(aRGB) == 2:
            global mode
            if not mode == "nothing":
                threadRunning = True
            else:
                threadRunning = False
            old_mode = mode
            mode = aRGB[0]
            setting = aRGB[1]
            if (mode == 'Rainbow'):
                global STEPS
                STEPS = float(setting)
                print "starting " + mode + "-service with " + str(STEPS) + " steps"
                if not mode == old_mode:
                        start_new_thread(self.rb.fader, (r,g,b))
            elif (mode == 'Music' or mode == 'Music1'):
                global MusicColor
                if mode == 'Music1':
                    MusicColor = setting.split('#')
                else:
                    MusicColor = setting
                print "starting " + mode + "-service with "  + str(MusicColor) + " color settings"
                if not mode == old_mode:
                        start_new_thread(self.ac.loop, ())
            elif (mode == "Flash"):
                global tempo
                tempo = float(setting)
                print "starting " + mode + "-service with a tempo of "  + str(tempo)
                if not mode == old_mode:
                        start_new_thread(self.ef.Flasher, ())
            elif (mode == "Strobe"):
                global tempo
                tempo = float(setting)
                print "starting " + mode + "-service with a tempo of "  + str(tempo)
                if not mode == old_mode:
                        start_new_thread(self.ef.Strober, ())

    def message_analyser(self,msg):
        RGBcolor = [0,0,0]
        if msg.startswith('#'):
            RGBcolor = self.lc.hex_to_rgb(msg)
        elif msg.startswith('rgb'):
            RGBcolor = [float(i) for i in msg[4:-1].split(',')]
        elif msg.startswith('hsl'):
            hslString = msg[4:-1].split(',')
            for i in range(3):
                if '%' in hslString[i]:
                    pos = hslString[i].index('%')
                    hslString[i] = hslString[i][:pos]
            RGBcolor = self.lc.hsl_to_rgb([float(i) for i in hslString])
        elif msg.count(',') == 1:
            RGBcolor = msg.split(',')
        elif msg.count(',') == 2 and not any(c.isalpha() for c in msg):
            RGBcolor = [float(i) for i in msg[4:-1].split(',')]
        else:
            RGBcolor = [0,0,0]
            print "Unsupported color model"
        return RGBcolor

    def on_close(self):
        # If client disconnects, remove him from the clients list
        self.clients.remove(self)

    def color_broadcaster(self,r,g,b):
        rgb = "rgb(" + str(r) + ", " + str(g) + ", " + str(b) + ")"
        self.send_message(rgb)

if __name__ == '__main__':
    if len(sys.argv) > 1:
        options['immediate_flush'] = False

    # 1. Create SockJSRouter
    BrokerRouter = sockjs.tornado.SockJSRouter(BrokerConnection, '/rgb')

    # 2. Create Tornado web.Application
    app = web.Application(BrokerRouter.urls)

    # 3. Make application listen on port
    app.listen(port)

    # 4. Every 1 second dump current client count
#    ioloop.PeriodicCallback(BrokerConnection.dump_stats, 1000).start()

    # 5. Start IOLoop
    ioloop.IOLoop.instance().start()

max_stars_repo_name: karaambaa/RGB-LED-Server

max_issues_repo_name: karaambaa/RGB-LED-Server

max_forks_repo_name: karaambaa/RGB-LED-Server

reason: import numpy,from scipy


text:
import time
import numpy as np
import pandas as pd
import os, sys
import pickle

from ..train.training import FSet, CV_fit, get_models, feat_comb_prepare
from ..train.training import reg_lib, clf_lib

class Benchmark:
    def __init__(self, mode='DDG_R', model_names=['RF_R'], Extra_model_dict=None, \
                 feat_combs=None, feat_alias=None, feats_lib=None, \
                 n_jobs=None, processed_dir=None, suffix=None, reloadTimestamp=None, \
                 FSet=FSet, clf_lib=clf_lib, reg_lib=reg_lib, feat_parallel=False, model_parallel=False):

        '''
            Parameters
            ----------
            mode: values of {'DG_R', 'DG_R', 'DDG_C'},
                    DG_R means regression of DG,
                    DDG_R means regression of DDG,
                    DDG_C means classification of DDG

            model_names:  list, eg: ['RF_R','XGB_R']

            Extra_model_dict: {k, model} like.

            feat_combs: list of feature combinations, if None, use all possible combinations.

            feat_alias: string or list of string, string should contains only the following letters:
                        'F', 'B', 'P', 's', 's', 'p'. Its meaning is：
                        {'F': 'mergefp', 'B': 'BoB', 'P': 'PhyChem_total',
                         's': 'SOAP', 'a': 'ACSF_local', 'p': 'PhyChem_local'}

            feats_lib: base feature type used to combinations, if feat_combs not None, feats_lib will be ignored.

            n_jobs: int or None, optional (default=None). The number of jobs to run in parallel.

            processed_dir: directory used to save generated files

            suffix: specified suffix of filename

            reloadTimestamp: used to reload saved files
        '''
        assert mode in ['DG_R', 'DDG_R',
                    'DDG_C'], 'mode should in DG_R/DDG_R/DDG_C'
        assert not (feat_alias and feat_combs), '\
            feat_combs and feat_alias are mutually exclusive and can \
            accept at most one of the parameters'

        self.mode = mode
        self.n_jobs = n_jobs
        self.processed_dir = processed_dir
        self.suffix = suffix
        self.reloadTimestamp = reloadTimestamp
        self.feat_parallel = feat_parallel
        self.model_parallel = model_parallel
        if not reloadTimestamp:
            self.LoadDataTime = time.strftime('%Y%m%d_%H%M%S', time.localtime(time.time()))
        else:
            self.LoadDataTime = reloadTimestamp

        if feat_combs:
            self.feat_combs, self.FSet = feat_combs, FSet
        else:
            self.feat_comb_prepare(feat_alias, feats_lib, FSet)
        self.get_models(model_names, Extra_model_dict, clf_lib, reg_lib)

    def feat_comb_prepare(self, feat_alias=None, feats_lib=None, FSet=FSet):
        '''
            Parameters
            ----------
            feat_combs: list of feature combinations, if None, use all possible combinations.

            feat_alias: string or list of string, string should contains only the following letters:
                        'F', 'B', 'P', 's', 's', 'p'. Its meaning is：
                        {'F': 'mergefp', 'B': 'BoB', 'P': 'PhyChem_total',
                         's': 'SOAP', 'a': 'ACSF_local', 'p': 'PhyChem_local'}

            feats_lib: base feature type used to combinations, if feat_combs not None, feats_lib will be ignored.
        '''
        global_feats_lib = ['mergefp', 'BoB', 'PhyChem_total', ] # , 'mergefp'， 'MACCSfp', 'Morganfp'
        local_feats_lib = ['SOAP', 'ACSF_local', 'PhyChem_local', ]

        self.feat_alias = feat_alias
        self.feats_lib = feats_lib
        self.feat_combs, self.FSet = feat_comb_prepare(feat_alias, feats_lib, FSet)

    def get_models(self, model_names=['RF_R'], Extra_model_dict=None, clf_lib=clf_lib, reg_lib=reg_lib):
        '''
            Parameters
            ----------
            model_names:  list, eg: ['RF_R','XGB_R']

            Extra_model_dict: {k, model} like.

            generate:
            -----------------
            self.model_lib: all model lib
            self.model_dict: used model
        '''
        self.Extra_model_dict = Extra_model_dict
        self.model_lib, self.model_dict = get_models(model_names,
                 self.mode, Extra_model_dict, clf_lib, reg_lib, n_jobs=self.n_jobs)
        self.model_names = list(self.model_dict.keys())

    def CV_fit(self, dataset, feat_combs=None, model_dict=None, inplace=False, **kwargs):
        '''
            Parameters
            ----------
            dataset:  dataset of torch_geometric.data.Data like, eg: ReactionDataset.data or SelectivityDataset.data

            feat_combs: list of feature combinations, if None, use all possible combinations.

            model_dict: {k, model} like.

            default values for kwargs:
                {'dataset': self.dataset, 'mode': self.mode,
                 'features_set': feat_combs, 'model_lib': self.model_lib,
                 'model_dict': model_dict, 'n_jobs': self.n_jobs, 'FSet': self.FSet}
        '''
        feat_combs = feat_combs if feat_combs else self.feat_combs
        model_dict = model_dict if model_dict else self.model_dict
        _kwargs = {'dataset': dataset, 'mode': self.mode,
                   'features_set': feat_combs, 'model_lib': self.model_lib,
                   'model_dict': model_dict, 'n_jobs': self.n_jobs, 'FSet': self.FSet,
                   'feat_parallel': self.feat_parallel, 'model_parallel': self.model_parallel}
        if kwargs:
            for k, v in kwargs:
                _kwargs[k] = v
        inplace = True if not (feat_combs or model_dict) else inplace

        model_res, y_pred_history, FSet, cv_result_history = CV_fit(**_kwargs)
        if inplace:
            self.summary = model_res
            self.y_pred_history = y_pred_history
            self.FSet = FSet
            self.CV_history = cv_result_history
        return model_res, y_pred_history, FSet, cv_result_history

    def save_to_pkl(self, obj, filename='Benchmark', storage_folder='Benchmark_pkl'):
        if os.path.isdir(r"%s/%s" % (self.processed_dir, storage_folder)) == False:
            os.makedirs(r"%s/%s" % (self.processed_dir, storage_folder))

        fn = r'%s/%s/%s_%s_%s.pkl' % (self.processed_dir, storage_folder,
                                             filename, self.suffix, self.LoadDataTime)
        print(fn)
        with open(fn,'wb') as f:
            pickle.dump(obj,f,0)

    def load_from_pkl(self, filename='Benchmark', storage_folder='Benchmark_pkl'):
        if os.path.isdir(r"%s/%s" % (self.processed_dir, storage_folder)) == False:
            os.makedirs(r"%s/%s" % (self.processed_dir, storage_folder))

        fn = r'%s/%s/%s_%s_%s.pkl' % (self.processed_dir, storage_folder,
                                             filename, self.suffix, self.reloadTimestamp)
        print(fn)
        with open(fn,'rb') as f:
            obj = pickle.load(f)
        return obj

    def __repr__(self):
        return '{}(n_feat_combs={}, n_models={})'.format(self.__class__.__name__,
                                        len(self.feat_combs), len(self.model_dict))

max_stars_repo_name: HongGroup/ChemSelML

max_issues_repo_name: HongGroup/ChemSelML

max_forks_repo_name: HongGroup/ChemSelML


text:
#!/usr/bin/python3
# number of output figures = 2

import matplotlib as mpl
import numpy as np
import scipy.interpolate

from helper.figure import Figure
import helper.plot
import helper.topo_opt

def main():
  h5Data = helper.topo_opt.readH5(
      "./data/topoOpt/results/650/thesis-2d-cantilever.h5")
  np.set_printoptions(threshold=np.nan)
  M1, M2 = 64, 32
  displacementFactor = 0.003

  displacement = np.reshape(h5Data["displacement"], (M2+1, M1+1, 2))
  displacement = np.transpose(displacement, (1, 0, 2))

  for q in range(2):
    fig = Figure.create(figsize=(5, 3), scale=0.6)
    ax = fig.gca()

    domainWidth, domainHeight = 2, 1
    macroCellWidth, macroCellHeight = domainWidth / M1, domainHeight / M2

    ax.add_patch(mpl.patches.Rectangle(
        (0, 0), domainWidth+0.016, domainHeight+0.016, edgecolor="none",
        facecolor=helper.plot.mixColors("mittelblau", 0.5),
        zorder=-10))

    microparams = h5Data["microparams"]["smart"]
    nn = np.array((M1, M2))
    VV = np.reshape(microparams[:,0] + microparams[:,1] -
                    microparams[:,0] * microparams[:,1], nn[::-1])
    XXYY = h5Data["nodes"][:,:2]
    VV = np.hstack((np.vstack((VV, VV[-1,:])),
                    np.reshape(np.append(VV[:,-1], VV[-1,-1]), (-1, 1))))
    VV[0,:]  = 0
    VV[-1,:] = 0
    VV[:,0]  = 0
    VV[:,-1] = 0
    VV = VV.flatten()
    triangulation = mpl.tri.Triangulation(XXYY[:,0], XXYY[:,1])
    ax.tricontourf(triangulation, VV, [0.25, 1], colors="hellblau")
    ax.tricontour(triangulation, VV, [0.25], colors="mittelblau")

    if q == 0:
      splines = [scipy.interpolate.RectBivariateSpline(
                    np.linspace(0, domainWidth,  M1+1),
                    np.linspace(0, domainHeight, M2+1),
                    displacement[:,:,t]) for t in range(2)]
      XXYYDisplaced = XXYY + displacementFactor * np.array(
          [splines[t](XXYY[:,0], XXYY[:,1], grid=False) for t in range(2)]).T
      triangulationDisplaced = mpl.tri.Triangulation(
          XXYYDisplaced[:,0], XXYYDisplaced[:,1])
      ax.tricontour(triangulationDisplaced, VV, [0.25],
                    colors=["k"], linestyles=["dashed"])

      ax.plot(domainWidth, 0, "k.", clip_on=False)
      helper.plot.plotArrow(ax, (domainWidth, 0), (domainWidth, -0.3))
      ax.text(domainWidth + 0.06, -0.15, r"$\force$", ha="left", va="center")

      ax.plot([0, 0], [0, domainHeight], "k-", clip_on=False)
      helper.plot.plotHatchedRectangle(
          ax, [-0.2, 0], [0.2, domainHeight], spacing=0.08, color="k")

      ax.text(0.6, 0.25, r"$\densglobal(\tilde{\*x}) = 1$",
              ha="center", va="center", rotation=20)
      ax.text(0.03, 0.45, r"$\densglobal(\tilde{\*x}) = 0$",
              ha="left", va="center")
      ax.text(domainWidth - 0.05, domainHeight - 0.05, r"$\objdomain$",
              ha="right", va="top")
    else:
      ax.plot([-0.2, -0.2], [0, domainHeight], "-", color=3*[0.99],
              clip_on=False)

    ax.set_aspect("equal")
    ax.set_xlim(-0.3, domainWidth)
    ax.set_ylim(-0.35, domainHeight)
    ax.set_axis_off()

    fig.save()

if __name__ == "__main__":
  main()

max_stars_repo_name: valentjn/defense

max_issues_repo_name: valentjn/defense

max_forks_repo_name: valentjn/defense

reason: import numpy,import scipy


text:
import numpy as np

class validate_raw_data:
    """
                            Class Name: validate_raw_data
                            Description: Validates the prediction data.
                            On Failure: Raise Exception

                            Written By: Murali Krishna Chintha
                            Version: 1.0
                            Revisions: None
    """

    def __init__(self, logger_obj, log_file):

        self.raw_prediction_validation_logs = log_file
        self.logger = logger_obj

    def column_validate(self, dataframe):
        """
                                Module Name: column_validate
                                Description: Validates the number of columns and the column names.
                                Input: dataframe
                                Output: True/False
                                On Failure: Raise Exception

                                Written By: Vaishnavi Ambati
                                Version: 1.0
                                Revisions: None
        """
        try:
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs, 'Entered column_validate module of validateRawPredictionData class.')
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs, 'Data Validation for columns has Stared.')

            columns = list(dataframe.columns)

            actual_columns = ['sku','national_inv','lead_time','in_transit_qty',
                              'forecast_3_month','forecast_6_month','forecast_9_month',
                              'sales_1_month','sales_3_month','sales_6_month',
                              'sales_9_month','min_bank','potential_issue',
                              'pieces_past_due','perf_6_month_avg','perf_12_month_avg',
                              'local_bo_qty','deck_risk','oe_constraint','ppap_risk',
                              'stop_auto_buy','rev_stop']

            # checking for column data types
            object_column_list = ['potential_issue',  'deck_risk','oe_constraint','ppap_risk',
                              'stop_auto_buy','rev_stop']

            float_column_list = ['national_inv','lead_time', 'sales_1_month', 'pieces_past_due',
                                 'perf_6_month_avg','local_bo_qty']

            object_columns = list(dataframe.dtypes[dataframe.dtypes == np.object].index)
            float_columns = list(dataframe.dtypes[dataframe.dtypes == np.float64].index)

            if object_columns == object_column_list and float_columns == float_column_list:
                data_type_validate = True
            else:
                data_type_validate = False

            if columns == actual_columns and data_type_validate:
                self.logger.write_log(self.raw_prediction_validation_logs, 'Columns in the data are valid.')
                self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'Exited column_validate module of validateRawPredictionData class.')

                return True
            elif columns == actual_columns and not data_type_validate:
                self.logger.write_log(self.raw_prediction_validation_logs, 'Columns in the dataset are valid but the column data types are not valid.')
                self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'Exited column_validate module of validateRawPredictionData class.')

                return False
            else:
                self.logger.write_log(self.raw_prediction_validation_logs, 'Columns in the data are not valid.')
                self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'Exited column_validate module of validateRawPredictionData class.')

                return False

        except Exception as e:
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'An Exception has occured in column_validate of validateRawPredicitonData class. The Exception is ' + str(e))
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs, 'Exting column_validate of validateRawPredicitonData class.')
            self.raw_prediction_validation_logs.to_csv("Logs\\Prediction Validation\\prediction_validation_logs.csv")

            raise Exception

    def precentage_column_null_value_check(self, dataframe):

        """
                                Module Name: entire_column_null_value_check
                                Description: Validates the number of columns and the column names.
                                Input: dataframe
                                Output: True/False
                                On Failure: Raise Exception

                                Written By: Vaishnavi Ambati
                                Version: 1.0
                                Revisions: None
        """

        self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs, 'Entered entire_column_null_value_check of validateRawPredicitonData class.')
        self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs, 'Null Values Check has started.')
        try:

            null_values_counts = dataframe.isnull().sum()
            length_of_dataframe = dataframe.shape[0]

            column_null_check = []

            for column, null_value_count in null_values_counts.items():
                if null_value_count >= 0.75 * length_of_dataframe:
                    self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,
                                                                                'The Column ' + str(column) + ' has only values. The Dataset Cannot be Accepted.')
                    self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,
                                                                                'Entire Row null value check has completed.')
                    self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,
                                                                                'Exting entireColumnNullValueCheck of validateRawPredicitonData class.')

                    column_null_check.append(True)
                else:
                    self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,
                                                                                'No column in the Dataset has complete null values.')
                    self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,
                                                                                'Exting entire_column_null_value_check of validateRawPredicitonData class.')

                    column_null_check.append(False)

            if False in column_null_check:
                return False
            else:
                return True

        except Exception as e:
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'An Exception has occured in entire_column_null_value_check of validateRawPredicitonData class. The Exception is ' + str(e))
            self.raw_prediction_validation_logs = self.logger.write_log(self.raw_prediction_validation_logs,'Exting entire_column_null_value_check of validateRawPredicitonData class.')
            self.raw_prediction_validation_logs.to_csv("Logs\\Prediction Validation\\prediction_validation_logs.csv")

            raise Exception