try: 
    import psyco
    psyco.full()
except:
    print "Impossible to import psyco library"
    pass
from importlibs import *



__author__ = "Andrea Maffezzoli"
__version__ = "$Revision: 0.35 - 10 March 2010$"
__name__ = "neurocorr"
__all__ = ['neurons_ac_binned', 'neurons_ac_isi', 'neurons_mi_binned', 'neurons_ccf']



def neurons_acf(neurons_ts_dict,
                bin = 0.1,
                delta = 1.0,
                start = 0.0,
                normalize = False,
                no_AC0 = False,
                smoothing_on = False,
                **smoothing_kwargs):
    """See also Neuroexplorer docs for the related explanations about
    how to calculate different kind of Auto-Correlation functions (ACF).
    Parameter normalize can be:
    - 'AC0', normalizing all ACFs values by dividing them to AC at zero delay;
    - 'all_ACs', normalizing all ACFs values by their sum;
    - 'bins_nr', normalizing all ACFs by the nr of bins (and so implicitly by
        the spike train's length);
    - False, if none normalization is required.
    - delta is the delay at which calculate the AC values.
    - start is the first delay at which give an AC value in output.
    - no_AC0 shows if to give in output thw AC value at zero delay
      (proportional to the variance).
    - smoothing is true if smoothing by convolution based on settings
      from file 'burst_detection_settings.txt' is required.
    
    Testing: to test ACFs compare their plots with the same from NeuroExplorer
    """

    acorr = dict.fromkeys(neurons_ts_dict)
    for neuron, ts_data in neurons_ts_dict.iteritems():
        if not len(ts_data): acorr[neuron] = numpy.asarray([]); continue
        acorr[neuron] = numpy.asarray(acf_OLD_VER(ts_data,
                                                  bin,
                                                  delta,
                                                  start = 0, 
                                                  two_side = False, 
                                                  no_bias = False,
                                                  no_AC_zero = no_AC0,
                                                  time_normalized = False)) # =False: normalization done below
        ###############################################################
        # Smooth AC values if requested, then convolve ACF with a pdf. 
        # An alternative: convolve original data, then calculate ACFs
        ###############################################################
        if smoothing_on: acorr[neuron] = convolve_with_pdf(acorr[neuron],
                                                           pdf = smoothing_kwargs['pdf'],
                                                           pdf_start = 0,
                                                           pdf_end = smoothing_kwargs['pdf_samples_nr'] * bin,
                                                           pdf_bin = bin)
        ################################################
        # Then we normalize AC values even if requested
        ################################################
        if normalize:
            if normalize   == 'all_ACs':  acorr_dividend = acorr[neuron].sum(dtype = numpy.float32)
            elif normalize == 'AC0':      acorr_dividend = acorr[neuron][0]
            elif normalize == 'duration': acorr_dividend = float(ts_data[-1] - ts_data[0])          
            if acorr_dividend: acorr[neuron] = numpy.divide(acorr[neuron], acorr_dividend)
    return acorr


def neurons_ccf(neurons_ts_dict,
                bin = 0.1,
                delta = 1.0,
                start = 0.0,
                normalize = False,
                smoothing_on = False,
                **smoothing_kwargs):
    """See also Neuroexplorer docs for the related explanations about
    how to calculate different kind of Auto-Correlation functions (ACF).
    Parameter normalize can be:
    - 'AC0', normalizing all ACFs values by dividing then to AC at zero delay;
    - 'all_ACs', normalizing all ACFs values by their sum;
    - 'bins_nr', normalizing all ACFs by the nr of bins (and so implicitly by
        the spike train's length);
    - 'False', if none normalization is required.
    Parameter delta is the delay at which calculate the AC values.
    Parameter start is the first delay at which give an AC value in output.
    Parameter smoothing is true if smoothing by convolution based on settings
    from file 'burst_detection_settings.txt' is required.
    """

    ccorr = {}
    for neuron1, neuron2 in set(combinations(neurons_ts_dict.iterkeys(), 2)):
        ccorr[(neuron1, neuron2)] = numpy.asarray(ccf(neurons_ts_dict[neuron1],
                                                      neurons_ts_dict[neuron2],
                                                      bin, delta,
                                                      start = 0,
                                                      no_bias = False,
                                                      normalized = False)) # False because normalization is done below
        
        ##############################
        # Normalizations of AC values
        ##############################
        if normalize:
            if normalize    == 'all_CCs':   ccorr_dividend = ccorr[(neuron1, neuron2)].sum(dtype = numpy.float32)
            elif normalize  == 'CC0':       ccorr_dividend = ccorr[(neuron1, neuron2)][0]
            elif normalize  == 'duration':  ccorr_dividend = float(neurons_ts_dict[neuron1][-1] - \
                                                                   neurons_ts_dict[neuron1][0])           
            if ccorr_dividend: ccorr[(neuron1, neuron2)] = numpy.divide(ccorr[(neuron1, neuron2)], ccorr_dividend)
        
        #####################################################################
        # Here we're convolving the ACF function with a given pdf.
        # An alternative would be to convolve at the begin the initial data,
        # whose the ACF is later calculated.
        #####################################################################
        if smoothing_on: ccorr[(neuron1, neuron2)] = convolve_with_pdf(ccorr[(neuron1, neuron2)],
                                                                       pdf = smoothing_kwargs['pdf'],
                                                                       pdf_start = 0,
                                                                       pdf_end = smoothing_kwargs['pdf_samples_nr'] * bin,
                                                                       pdf_bin = bin)
        
    return ccorr



def neurons_ac_isi(isi_series, delta = 0.0, normalize = True):
    """Calculating Auto-Correlation function values on ISI data"""

    acsi = {}
    for neuron, single_isi in isi_series.items():
        acsi_temp = ac_isi(single_isi, delta)

        ########################################################################
        # To normalize ACs at all delays to zero delay AC, divide them to AC[0]
        ########################################################################
        if normalize: acsi[neuron] = numpy.divide(acsi_temp, float(acsi_temp[0]))
        else: acsi[neuron] = numpy.asarray(acsi_temp)
    
    return acsi


def neurons_mi_binned(histogram_countsperbin, bin = 1.0, delta = 0.0, normalize = True):
    """See also the related explanations for Auto-Correlation functions"""

    misb = {}
    for neuron, histo in histogram_countsperbin.items():
        if any(histo):
            ######################################################
            # In the following using as input of AC function:
            # - the spikes counts per bin by the total spikes nr.
            ######################################################
            histogram_data = histo

            #############################################################################
            # - the discrete distribution of probability of spikes per bin,
            #   calculated by dividing the spikes counts per bin by the total spikes nr.
            #############################################################################
            #histogram_data = numpy.divide(histogram_countsperbin, float(len(histo)))
            
            ######################################################################
            # - the distribution of spikes per second, calculated by dividing the 
            #   spikes counts per bin by the total spikes nr per bin lenght.
            ######################################################################
            #histogram_data = numpy.divide(histo, numpy.sum(histo) * bin)
            
            misb_temp = mi_binned(histogram_data, bin, delta)
            
            ##############################################################################################
            # To normalize MI values for all delays with MI value at zero delay, divide the MIs for MI[0]
            ##############################################################################################
            if normalize: misb[neuron] = numpy.divide(misb_temp, float(misb_temp[0]))
            else: misb[neuron] = numpy.asarray(misb_temp)

    return misb
