""" 

Copyright (c) 2016, Jose Dolz .All rights reserved.

Redistribution and use in source and binary forms, with or without modification,

are permitted provided that the following conditions are met:

    1. Redistributions of source code must retain the above copyright notice,

       this list of conditions and the following disclaimer.

    2. Redistributions in binary form must reproduce the above copyright notice,

       this list of conditions and the following disclaimer in the documentation

       and/or other materials provided with the distribution.

    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,

    EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES

    OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND

    NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT

    HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,

    WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING

    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR

    OTHER DEALINGS IN THE SOFTWARE.

NOTES: There are still some functionalities to be implemented.

    - Add pooling layer in 3D

    - Add more activation functions

    - Add more optimizers (ex. Adam)

Jose Dolz. Dec, 2016.

email: jose.dolz.upv@gmail.com

LIVIA Department, ETS, Montreal.

"""

import numpy

import numpy as np

import theano

import theano.tensor as T

from theano.tensor.nnet import conv

import random

from math import floor

from math import ceil

from Modules.General.Utils import computeReceptiveField

from Modules.General.Utils import extendLearningRateToParams

from Modules.General.Utils import extractCenterFeatMaps

from Modules.General.Utils import getCentralVoxels

from Modules.General.Utils import getWeightsSet

import LiviaNet3DConvLayer

import LiviaSoftmax

import pdb

#####################################################

# ------------------------------------------------- #

##  ##  ##  ##  ##   LIVIANET 3D   ##  ##  ##  ##  ##

# ------------------------------------------------- #

#####################################################

class LiviaNet3D(object):

    def __init__(self):

        # --- containers for Theano compiled functions ----

        self.networkModel_Train = ""

        self.networkModel_Test = ""

        # --- shared variables will be stored in the following variables ----

        self.trainingData_x = ""

        self.testingData_x = ""

        self.trainingData_y = ""

        self.lastLayer = ""

        self.networkLayers = []

        self.intermediate_ConnectedLayers = []

        self.networkName = ""

        self.folderName = ""

        self.cnnLayers = []

        self.n_classes = -1

        self.sampleSize_Train = []

        self.sampleSize_Test = []

        self.kernel_Shapes = []

        self.pooling_scales = []

        self.dropout_Rates = []

        self.activationType = -1

        self.weight_Initialization = -1

        self.dropoutRates = []

        self.batch_Size = -1

        self.receptiveField = 0

        self.initialLearningRate = "" 

        self.learning_rate = theano.shared(np.cast["float32"](0.01))

        # Symbolic variables,

        self.inputNetwork_Train = None  

        self.inputNetwork_Test = None

        self.L1_reg_C = 0

        self.L2_reg_C = 0

        self.costFunction = 0

        # Params for optimizers

        self.initialMomentum = "" 

        self.momentum = theano.shared(np.cast["float32"](0.))

        self.momentumNormalized = 0

        self.momentumType = 0

        self.vel_Momentum = [] 

        self.rho_RMSProp = 0

        self.epsilon_RMSProp = 0

        self.params_RmsProp = [] 

        self.numberOfEpochsTrained = 0

        self.applyBatchNorm = ""

        self.numberEpochToApplyBatchNorm = 0

        self.softmax_Temp = 1.0

        self.centralVoxelsTrain = ""

        self.centralVoxelsTest = ""

    # -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    """ ####### Function to generate the network architecture  ######### """

    def generateNetworkLayers(self,

                            cnnLayers,

                            kernel_Shapes,

                            maxPooling_Layer,

                            sampleShape_Train,

                            sampleShape_Test,

                            inputSample_Train,

                            inputSample_Test,

                            layersToConnect):

        rng = np.random.RandomState(24575)

        # Define inputs for first layers (which will be re-used for next layers)

        inputSampleToNextLayer_Train = inputSample_Train

        inputSampleToNextLayer_Test = inputSample_Test

        inputSampleToNextLayerShape_Train = sampleShape_Train

        inputSampleToNextLayerShape_Test = sampleShape_Test

        # Get the convolutional layers

        numLayers = len(kernel_Shapes)

        numberCNNLayers = []

        numberFCLayers = []

        for l_i in range(1,len(kernel_Shapes)):

            if len(kernel_Shapes[l_i]) == 3:

                numberCNNLayers = l_i + 1

        numberFCLayers = numLayers - numberCNNLayers

        ######### -------------- Generate the convolutional layers --------------   #########

        # Some checks

        if self.weight_Initialization_CNN == 2:

            if len(self.weightsTrainedIdx) <> numberCNNLayers:

                print(" ... WARNING!!!! Number of indexes specified for trained layers does not correspond with number of conv layers in the created architecture...")

        if self.weight_Initialization_CNN == 2:

            weightsNames = getWeightsSet(self.weightsFolderName, self.weightsTrainedIdx)

        for l_i in xrange(0, numberCNNLayers) :

            # Get properties of this layer

            # The second element is the number of feature maps of previous layer

            currentLayerKernelShape = [cnnLayers[l_i], inputSampleToNextLayerShape_Train[1]] +  kernel_Shapes[l_i] 

            # If weights are going to be initialized from other pre-trained network they should be loaded in this stage

            # Otherwise

            weights = []

            if self.weight_Initialization_CNN == 2:

                weights = np.load(weightsNames[l_i])

            maxPoolingParameters = []

            dropoutRate = 0.0

            myLiviaNet3DConvLayer = LiviaNet3DConvLayer.LiviaNet3DConvLayer(rng,

                                                                            l_i,

                                                                            inputSampleToNextLayer_Train,

                                                                            inputSampleToNextLayer_Test,

                                                                            inputSampleToNextLayerShape_Train,

                                                                            inputSampleToNextLayerShape_Test,

                                                                            currentLayerKernelShape,

                                                                            self.applyBatchNorm,

                                                                            self.numberEpochToApplyBatchNorm,

                                                                            maxPoolingParameters,

                                                                            self.weight_Initialization_CNN,

                                                                            weights,

                                                                            self.activationType,

                                                                            dropoutRate

                                                                            )

            self.networkLayers.append(myLiviaNet3DConvLayer)

            # Just for printing

            inputSampleToNextLayer_Train_Old = inputSampleToNextLayerShape_Train

            inputSampleToNextLayer_Test_Old  = inputSampleToNextLayerShape_Test

            # Update inputs for next layer

            inputSampleToNextLayer_Train = myLiviaNet3DConvLayer.outputTrain

            inputSampleToNextLayer_Test  = myLiviaNet3DConvLayer.outputTest    

            inputSampleToNextLayerShape_Train = myLiviaNet3DConvLayer.outputShapeTrain

            inputSampleToNextLayerShape_Test  = myLiviaNet3DConvLayer.outputShapeTest

            print(" ----- (Training) Input shape: {}  ---> Output shape: {}  ||  kernel shape {}".format(inputSampleToNextLayer_Train_Old,inputSampleToNextLayerShape_Train, currentLayerKernelShape))

            print(" ----- (Testing) Input shape: {}   ---> Output shape: {}".format(inputSampleToNextLayer_Test_Old,inputSampleToNextLayerShape_Test))

        ######### -------------- Create the intermediate (i.e. multi-scale) connections from conv layers to FCN ----------------- ##################

        featMapsInFullyCN = inputSampleToNextLayerShape_Train[1]

        [featMapsInFullyCN, 

        inputToFullyCN_Train,

        inputToFullyCN_Test] = self.connectIntermediateLayers(layersToConnect,

                                                              inputSampleToNextLayer_Train,

                                                              inputSampleToNextLayer_Test,

                                                              featMapsInFullyCN)

        ######### --------------  Generate the Fully Connected Layers  ----------------- ##################

        # Define inputs

        inputFullyCNShape_Train = [self.batch_Size, featMapsInFullyCN] + inputSampleToNextLayerShape_Train[2:5]

        inputFullyCNShape_Test = [self.batch_Size, featMapsInFullyCN] + inputSampleToNextLayerShape_Test[2:5]

        # Kamnitsas applied padding and mirroring to the images when kernels in FC layers were larger than 1x1x1.

        # For this current work, we employed kernels of this size (i.e. 1x1x1), so there is no need to apply padding or mirroring.

        # TODO. Check

        print(" --- Starting to create the fully connected layers....")

        for l_i in xrange(numberCNNLayers, numLayers) :

            numberOfKernels = cnnLayers[l_i]

            kernel_shape = [kernel_Shapes[l_i][0],kernel_Shapes[l_i][0],kernel_Shapes[l_i][0]]

            currentLayerKernelShape = [cnnLayers[l_i], inputFullyCNShape_Train[1]] +  kernel_shape

            # If weights are going to be initialized from other pre-trained network they should be loaded in this stage

            # Otherwise

            weights = []

            applyBatchNorm = True

            epochsToApplyBatchNorm = 60

            maxPoolingParameters = []

            dropoutRate = self.dropout_Rates[l_i-numberCNNLayers]

            myLiviaNet3DFullyConnectedLayer = LiviaNet3DConvLayer.LiviaNet3DConvLayer(rng,

                                                                            l_i,

                                                                            inputToFullyCN_Train,

                                                                            inputToFullyCN_Test,

                                                                            inputFullyCNShape_Train,

                                                                            inputFullyCNShape_Test,

                                                                            currentLayerKernelShape,

                                                                            self.applyBatchNorm,

                                                                            self.numberEpochToApplyBatchNorm,

                                                                            maxPoolingParameters,

                                                                            self.weight_Initialization_FCN,

                                                                            weights,

                                                                            self.activationType,

                                                                            dropoutRate

                                                                            )

            self.networkLayers.append(myLiviaNet3DFullyConnectedLayer)

            # Just for printing

            inputFullyCNShape_Train_Old = inputFullyCNShape_Train

            inputFullyCNShape_Test_Old  = inputFullyCNShape_Test

            # Update inputs for next layer

            inputToFullyCN_Train = myLiviaNet3DFullyConnectedLayer.outputTrain

            inputToFullyCN_Test = myLiviaNet3DFullyConnectedLayer.outputTest    

            inputFullyCNShape_Train = myLiviaNet3DFullyConnectedLayer.outputShapeTrain

            inputFullyCNShape_Test = myLiviaNet3DFullyConnectedLayer.outputShapeTest

            # Print

            print(" ----- (Training) Input shape: {}  ---> Output shape: {}  ||  kernel shape {}".format(inputFullyCNShape_Train_Old,inputFullyCNShape_Train, currentLayerKernelShape))

            print(" ----- (Testing) Input shape: {}   ---> Output shape: {}".format(inputFullyCNShape_Test_Old,inputFullyCNShape_Test))

        ######### -------------- Do Classification layer  ----------------- ##################

        # Define kernel shape for classification layer

        featMaps_LastLayer = self.cnnLayers[-1]

        filterShape_ClassificationLayer = [self.n_classes, featMaps_LastLayer, 1, 1, 1]

        # Define inputs and shapes for the classification layer

        inputImageClassificationLayer_Train = inputToFullyCN_Train

        inputImageClassificationLayer_Test = inputToFullyCN_Test

        inputImageClassificationLayerShape_Train = inputFullyCNShape_Train

        inputImageClassificationLayerShape_Test = inputFullyCNShape_Test

        print(" ----- (Classification layer) kernel shape {}".format(filterShape_ClassificationLayer))

        classification_layer_Index = l_i

        weights = []

        applyBatchNorm = True

        epochsToApplyBatchNorm = 60

        maxPoolingParameters = []

        dropoutRate = self.dropout_Rates[len(self.dropout_Rates)-1]

        softmaxTemperature = 1.0

        myLiviaNet_ClassificationLayer = LiviaSoftmax.LiviaSoftmax(rng,

                                                                   classification_layer_Index,

                                                                   inputImageClassificationLayer_Train,

                                                                   inputImageClassificationLayer_Test,

                                                                   inputImageClassificationLayerShape_Train,

                                                                   inputImageClassificationLayerShape_Test,

                                                                   filterShape_ClassificationLayer,

                                                                   self.applyBatchNorm,

                                                                   self.numberEpochToApplyBatchNorm,

                                                                   maxPoolingParameters,

                                                                   self.weight_Initialization_FCN,

                                                                   weights,

                                                                   0, #self.activationType,

                                                                   dropoutRate,

                                                                   softmaxTemperature

                                                                   )

        self.networkLayers.append(myLiviaNet_ClassificationLayer)

        self.lastLayer = myLiviaNet_ClassificationLayer

        print(" ----- (Training) Input shape: {}  ---> Output shape: {}  ||  kernel shape {}".format(inputImageClassificationLayerShape_Train,myLiviaNet_ClassificationLayer.outputShapeTrain, filterShape_ClassificationLayer))

        print(" ----- (Testing) Input shape:  {}  ---> Output shape: {}".format(inputImageClassificationLayerShape_Test,myLiviaNet_ClassificationLayer.outputShapeTest))

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    def updateLayersMatricesBatchNorm(self):

        for l_i in xrange(0, len(self.networkLayers) ) :

            self.networkLayers[l_i].updateLayerMatricesBatchNorm()

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    """ Function that connects intermediate layers to the input of the first fully connected layer 

        This is done for multi-scale features """    

    def connectIntermediateLayers(self,

                                  layersToConnect,

                                  inputSampleInFullyCN_Train,

                                  inputSampleInFullyCN_Test,

                                  featMapsInFullyCN):

        centralVoxelsTrain = self.centralVoxelsTrain

        centralVoxelsTest = self.centralVoxelsTest

        for l_i in layersToConnect :

            currentLayer = self.networkLayers[l_i]

            output_train = currentLayer.outputTrain

            output_trainShape = currentLayer.outputShapeTrain

            output_test = currentLayer.outputTest

            output_testShape = currentLayer.outputShapeTest

            # Get the middle part of feature maps at intermediate levels to make them of the same shape at the beginning of the

            # first fully connected layer

            featMapsCenter_Train = extractCenterFeatMaps(output_train, output_trainShape, centralVoxelsTrain)

            featMapsCenter_Test  = extractCenterFeatMaps(output_test, output_testShape, centralVoxelsTest)

            featMapsInFullyCN = featMapsInFullyCN + currentLayer._numberOfFeatureMaps

            inputSampleInFullyCN_Train = T.concatenate([inputSampleInFullyCN_Train, featMapsCenter_Train], axis=1)

            inputSampleInFullyCN_Test = T.concatenate([inputSampleInFullyCN_Test, featMapsCenter_Test], axis=1)

        return [featMapsInFullyCN, inputSampleInFullyCN_Train, inputSampleInFullyCN_Test]

    #############   Functions for OPTIMIZERS ################# 

    def getUpdatesOfTrainableParameters(self, cost, paramsTraining, numberParamsPerLayer) :

        # Optimizers

        def SGD():

            print (" --- Optimizer: Stochastic gradient descent (SGD)")

            updates = self.updateParams_SGD(cost, paramsTraining, numberParamsPerLayer)

            return updates

        def RMSProp():

            print (" --- Optimizer: RMS Prop")

            updates = self.updateParams_RMSProp(cost, paramsTraining, numberParamsPerLayer)

            return updates

        # TODO. Include more optimizers here

        optionsOptimizer = {0 : SGD,

                            1 : RMSProp}

        updates = optionsOptimizer[self.optimizerType]()

        return updates

    """ # Optimizers:

    # More optimizers in : https://github.com/Lasagne/Lasagne/blob/master/lasagne/updates.py """

    # ========= Update the trainable parameters using Stocastic Gradient Descent ===============

    def updateParams_SGD(self, cost, paramsTraining, numberParamsPerLayer) :

        # Create a list of gradients for all model parameters

        grads = T.grad(cost, paramsTraining)

        # Get learning rates for each param

        #learning_rates = extendLearningRateToParams(numberParamsPerLayer,self.learning_rate)

        self.vel_Momentum = []

        updates = []

        constantForCurrentGradientUpdate = 1.0 - self.momentum*self.momentumNormalized 

        #for param, grad, lrate  in zip(paramsTraining, grads, learning_rates) :

        for param, grad  in zip(paramsTraining, grads) :

            v = theano.shared(param.get_value()*0., broadcastable=param.broadcastable)

            self.vel_Momentum.append(v)

            stepToGradientDirection = constantForCurrentGradientUpdate*self.learning_rate*grad

            newVel = self.momentum * v - stepToGradientDirection

            if self.momentumType == 0 : 

                updateToParam = newVel

            else : 

                updateToParam = self.momentum*newVel - stepToGradientDirection

            updates.append((v, newVel)) 

            updates.append((param, param + updateToParam))

        return updates

    # ========= Update the trainable parameters using RMSProp ===============

    def updateParams_RMSProp(self, cost, paramsTraining, numberParamsPerLayer) : 

        # Original code: https://gist.github.com/Newmu/acb738767acb4788bac3

        # epsilon=1e-4 in paper.

        # Kamnitsas reported NaN values in cost function when employing this value.

        # Worked ok with epsilon=1e-6.

        grads = T.grad(cost, paramsTraining)

        # Get learning rates for each param

        #learning_rates = extendLearningRateToParams(numberParamsPerLayer,self.learning_rate)

        self.params_RmsProp = []

        self.vel_Momentum = []

        updates = []

        constantForCurrentGradientUpdate = 1.0 - self.momentum*self.momentumNormalized 

        # Using theano constant to prevent upcasting of float32

        one = T.constant(1)

        for param, grad in zip(paramsTraining, grads):

            accu = theano.shared(param.get_value()*0., broadcastable=param.broadcastable)

            self.params_RmsProp.append(accu) 

            v = theano.shared(param.get_value()*0., broadcastable=param.broadcastable) 

            self.vel_Momentum.append(v)

            accu_new = self.rho_RMSProp * accu + (one - self.rho_RMSProp) * T.sqr(grad)

            numGradStep = self.learning_rate * grad

            denGradStep = T.sqrt(accu_new + self.epsilon_RMSProp)

            stepToGradientDirection = constantForCurrentGradientUpdate*(numGradStep /denGradStep) 

            newVel = self.momentum * v - stepToGradientDirection

            if self.momentumType == 0 : 

                updateToParam = newVel

            else : 

                updateToParam = self.momentum*newVel - stepToGradientDirection

            updates.append((accu, accu_new))

            updates.append((v, newVel)) 

            updates.append((param, param + updateToParam))

        return updates

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    """ ------ Get trainable parameters --------- """

    def getTrainable_Params(_self):

        trainable_Params = []

        numberTrain_ParamsLayer = [] 

        for l_i in xrange(0, len(_self.networkLayers) ) :

            trainable_Params = trainable_Params + _self.networkLayers[l_i].params

            numberTrain_ParamsLayer.append(_self.networkLayers[l_i].numberOfTrainableParams) # TODO: Get this directly as len(_self.networkLayers[l_i].params)

        return trainable_Params,numberTrain_ParamsLayer

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    def initTrainingParameters(self,

                               costFunction,

                               L1_reg_C,

                               L2_reg_C,

                               learning_rate,

                               momentumType,

                               momentumValue,

                               momentumNormalized,

                               optimizerType,

                               rho_RMSProp,

                               epsilon_RMSProp

                               ) :

        print(" ------- Initializing network training parameters...........")

        self.numberOfEpochsTrained = 0

        self.L1_reg_C = L1_reg_C

        self.L2_reg_C = L2_reg_C

        # Set Learning rate and store the last epoch where it was modified

        self.initialLearningRate = learning_rate

        # TODO: Check the shared variables from learning rates

        self.learning_rate.set_value(self.initialLearningRate[0])          

        # Set momentum type and values

        self.momentumType = momentumType

        self.initialMomentumValue = momentumValue

        self.momentumNormalized = momentumNormalized

        self.momentum.set_value(self.initialMomentumValue)

        # Optimizers

        if (optimizerType == 2):

            optimizerType = 1

        def SGD():

            print (" --- Optimizer: Stochastic gradient descent (SGD)")

            self.optimizerType = optimizerType

        def RMSProp():

            print (" --- Optimizer: RMS Prop")

            self.optimizerType = optimizerType

            self.rho_RMSProp = rho_RMSProp

            self.epsilon_RMSProp = epsilon_RMSProp

        # TODO. Include more optimizers here

        optionsOptimizer = {0 : SGD,

                            1 : RMSProp}

        optionsOptimizer[optimizerType]()

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

    def updateParams_BatchNorm(self) : 

        updatesForBnRollingAverage = []

        for l_i in xrange(0, len(self.networkLayers) ) :

            currentLayer = self.networkLayers[l_i]

            updatesForBnRollingAverage.extend( currentLayer.getUpdatesForBnRollingAverage() ) 

        return updatesForBnRollingAverage

    # ------------------------------------------------------------------------------------ #

    # ---------------------------     Compile the Theano functions     ------------------- #

    # ------------------------------------------------------------------------------------ #

    def compileTheanoFunctions(self):

        print(" ----------------- Starting compilation process ----------------- ")        

        # ------- Create and initialize sharedVariables needed to compile the training function ------ #

        # -------------------------------------------------------------------------------------------- #

        # For training 

        self.trainingData_x = theano.shared(np.zeros([1,1,1,1,1], dtype="float32"), borrow = True)

        self.trainingData_y = theano.shared(np.zeros([1,1,1,1], dtype="float32") , borrow = True)  

        # For testing 

        self.testingData_x = theano.shared(np.zeros([1,1,1,1,1], dtype="float32"), borrow = True)

        x_Train = self.inputNetwork_Train

        x_Test  = self.inputNetwork_Test

        y_Train = T.itensor4('y')

        # Allocate symbolic variables for the data

        index_Train = T.lscalar()

        index_Test  = T.lscalar()

        # ------- Needed to compile the training function ------ #

        # ------------------------------------------------------ #

        trainingData_y_CastedToInt   = T.cast( self.trainingData_y, 'int32') 

        # To accomodate the weights in the cost function to account for class imbalance

        weightsOfClassesInCostFunction = T.fvector()  

        weightPerClass = T.fvector() 

        # --------- Get trainable parameters (to be fit by gradient descent) ------- #

        # -------------------------------------------------------------------------- #

        [paramsTraining, numberParamsPerLayer] = self.getTrainable_Params()

        # ------------------ Define the cost function --------------------- #

        # ----------------------------------------------------------------- #

        def negLogLikelihood():

            print (" --- Cost function: negativeLogLikelihood")

            costInLastLayer = self.lastLayer.negativeLogLikelihoodWeighted(y_Train,weightPerClass)

            return costInLastLayer

        def NotDefined():

            print (" --- Cost function: Not defined!!!!!! WARNING!!!")

        optionsCostFunction = {0 : negLogLikelihood,

                               1 : NotDefined}

        costInLastLayer = optionsCostFunction[self.costFunction]()

        # --------------------------- Get costs --------------------------- #

        # ----------------------------------------------------------------- #

        # Get L1 and L2 weights regularization

        costL1 = 0

        costL2 = 0

        # Compute the costs

        for l_i in xrange(0, len(self.networkLayers)) :    

                costL1 += abs(self.networkLayers[l_i].W).sum()

                costL2 += (self.networkLayers[l_i].W ** 2).sum()

        # Add also the cost of the last layer                     

        cost = (costInLastLayer

                + self.L1_reg_C * costL1

                + self.L2_reg_C * costL2)

        # --------------------- Include all trainable parameters in updates (for optimization) ---------------------- #

        # ----------------------------------------------------------------------------------------------------------- #

        updates = self.getUpdatesOfTrainableParameters(cost, paramsTraining, numberParamsPerLayer)

        # --------------------- Include batch normalization params ---------------------- #

        # ------------------------------------------------------------------------------- #

        updates = updates + self.updateParams_BatchNorm()

        # For the testing function we need to get the Feature maps activations

        featMapsActivations = []

        lower_act = 0

        upper_act = 9999

        # TODO: Change to output_Test

        for l_i in xrange(0,len(self.networkLayers)):

            featMapsActivations.append(self.networkLayers[l_i].outputTest[:, lower_act : upper_act, :, :, :])

        # For the last layer get the predicted probabilities (p_y_given_x_test)

        featMapsActivations.append(self.lastLayer.p_y_given_x_test)

        # --------------------- Preparing data to compile the functions ---------------------- #

        # ------------------------------------------------------------------------------------ #

        givensDataSet_Train = { x_Train: self.trainingData_x[index_Train * self.batch_Size: (index_Train + 1) * self.batch_Size],

                                y_Train: trainingData_y_CastedToInt[index_Train * self.batch_Size: (index_Train + 1) * self.batch_Size],

                                weightPerClass: weightsOfClassesInCostFunction }

        givensDataSet_Test  = { x_Test: self.testingData_x[index_Test * self.batch_Size: (index_Test + 1) * self.batch_Size] }

        print(" ...Compiling the training function...")

        self.networkModel_Train = theano.function(

                                    [index_Train, weightsOfClassesInCostFunction],

                                    #[cost] + self.lastLayer.doEvaluation(y_Train),

                                    [cost],

                                    updates=updates,

                                    givens = givensDataSet_Train

                                    )

        print(" ...The training function was compiled...")

        #self.getProbabilities = theano.function(

                         #[index],

                         #self.lastLayer.p_y_given_x_Train,

                         #givens={

                            #x: self.trainingData_x[index * _self.batch_size: (index + 1) * _self.batch_size]

                         #}

         #)

        print(" ...Compiling the testing function...")

        self.networkModel_Test = theano.function(

                                  [index_Test],

                                  featMapsActivations,

                                  givens = givensDataSet_Test

                                  )

        print(" ...The testing function was compiled...")

# -------------------------------------------------------------------- END Function ------------------------------------------------------------------- #

####### Function to generate the CNN #########

    def createNetwork(self,

                      networkName, 

                      folderName,

                      cnnLayers,

                      kernel_Shapes,

                      intermediate_ConnectedLayers,

                      n_classes,

                      sampleSize_Train,

                      sampleSize_Test,

                      batch_Size,

                      applyBatchNorm,

                      numberEpochToApplyBatchNorm,

                      activationType,

                      dropout_Rates,

                      pooling_Params,

                      weights_Initialization_CNN,

                      weights_Initialization_FCN,

                      weightsFolderName,

                      weightsTrainedIdx,

                      softmax_Temp

                      ):

        # ============= Model Parameters Passed as arguments ================

        # Assign parameters:

        self.networkName = networkName

        self.folderName = folderName

        self.cnnLayers = cnnLayers

        self.n_classes = n_classes

        self.kernel_Shapes = kernel_Shapes

        self.intermediate_ConnectedLayers = intermediate_ConnectedLayers

        self.pooling_scales = pooling_Params

        self.dropout_Rates = dropout_Rates

        self.activationType = activationType

        self.weight_Initialization_CNN = weights_Initialization_CNN

        self.weight_Initialization_FCN = weights_Initialization_FCN

        self.weightsFolderName = weightsFolderName

        self.weightsTrainedIdx = weightsTrainedIdx

        self.batch_Size = batch_Size

        self.sampleSize_Train = sampleSize_Train

        self.sampleSize_Test = sampleSize_Test

        self.applyBatchNorm = applyBatchNorm

        self.numberEpochToApplyBatchNorm = numberEpochToApplyBatchNorm

        self.softmax_Temp = softmax_Temp

        # Compute the CNN receptive field

        stride = 1;

        self.receptiveField = computeReceptiveField(self.kernel_Shapes, stride)

        # --- Size of Image samples ---

        self.sampleSize_Train = sampleSize_Train

        self.sampleSize_Test = sampleSize_Test

        ## --- Batch Size ---

        self.batch_Size = batch_Size

        # ======== Calculated Attributes =========

        self.centralVoxelsTrain = getCentralVoxels(self.sampleSize_Train, self.receptiveField) 

        self.centralVoxelsTest = getCentralVoxels(self.sampleSize_Test, self.receptiveField) 

        #==============================

        rng = numpy.random.RandomState(23455)

        # Transfer to LIVIA NET

        self.sampleSize_Train = sampleSize_Train

        self.sampleSize_Test = sampleSize_Test

        # --------- Now we build the model -------- #

        print("...[STATUS]: Building the Network model...")

        # Define the symbolic variables used as input of the CNN

        # start-snippet-1

        # Define tensor5

        tensor5 = T.TensorType(dtype='float32', broadcastable=(False, False, False, False, False))

        self.inputNetwork_Train = tensor5() 

        self.inputNetwork_Test = tensor5()

        # Define input shapes to the netwrok

        inputSampleShape_Train = (self.batch_Size, 1, self.sampleSize_Train[0], self.sampleSize_Train[1], self.sampleSize_Train[2])

        inputSampleShape_Test = (self.batch_Size, 1, self.sampleSize_Test[0], self.sampleSize_Test[1], self.sampleSize_Test[2])

        print (" - Shape of input subvolume (Training): {}".format(inputSampleShape_Train))

        print (" - Shape of input subvolume (Testing): {}".format(inputSampleShape_Test))

        inputSample_Train = self.inputNetwork_Train

        inputSample_Test = self.inputNetwork_Test

        # TODO change cnnLayers name by networkLayers

        self.generateNetworkLayers(cnnLayers,

                                   kernel_Shapes,

                                   self.pooling_scales,

                                   inputSampleShape_Train,

                                   inputSampleShape_Test,

                                   inputSample_Train,

                                   inputSample_Test,

                                   intermediate_ConnectedLayers)      

    # Release Data from GPU

    def releaseGPUData(self) :

        # GPU NOTE: Remove the input values to avoid copying data to the GPU

        # Image Data

        self.trainingData_x.set_value(np.zeros([1,1,1,1,1], dtype="float32"))

        self.testingData_x.set_value(np.zeros([1,1,1,1,1], dtype="float32"))

        # Labels

        self.trainingData_y.set_value(np.zeros([1,1,1,1], dtype="float32"))