Python scipy.linspace() Examples

def sim_and_plot(G, tau, gamma, rho, tmax, tcount, ax):
    t, S, I = EoN.fast_SIS(G, tau, gamma, rho = rho, tmax = tmax)
    report_times = scipy.linspace(0, tmax, tcount)
    I = EoN.subsample(report_times, t, I)
    ax.plot(report_times, I/N, color='grey', linewidth=5, alpha=0.3)
    
    t, S, I, = EoN.SIS_heterogeneous_meanfield_from_graph(G, tau, gamma, rho=rho, 
                                                    tmax=tmax, tcount=tcount)
    ax.plot(t, I/N, '--')    
    t, S, I = EoN.SIS_compact_pairwise_from_graph(G, tau, gamma, rho=rho,
                                                    tmax=tmax, tcount=tcount)
    ax.plot(t, I/N)
 
    t, S, I = EoN.SIS_homogeneous_pairwise_from_graph(G, tau, gamma, rho=rho, 
                                                    tmax=tmax, tcount=tcount)
    ax.plot(t, I/N, '-.') 

def SIS_process(G, degree_prob, tmax, tau, gamma):
    N=G.order() 
    plt.figure(5)
    plt.clf()
    plt.figure(6)
    plt.clf()
    for index, starting_node in enumerate([x*N/10. for x in range(10)]):
        plt.figure(5)
        t, S, I = EoN.fast_SIS(G, tau, gamma, initial_infecteds = [starting_node], tmax = tmax)
        #print(I[-1])
        subt = scipy.linspace(0, tmax, 501)
        subI = EoN.subsample(subt, t, I)
        plt.plot(subt,subI)
        if I[-1]>100:
            plt.figure(6)
            shift = EoN.get_time_shift(t, I, 1000)
            plt.plot(subt-shift, subI)
    plt.figure(5)
    plt.savefig('sw_SIS_epi_N{}_p{}_k{}_tau{}.pdf'.format(N, p, k, tau))
    plt.figure(6)
    plt.savefig('sw_SIS_epi_N{}_p{}_k{}_tau{}_shifted.pdf'.format(N, p, k, tau)) 

def m_circles(mags, phase_min=-359.75, phase_max=-0.25):
    """Constant-magnitude contours of the function Gcl = Gol/(1+Gol), where
    Gol is an open-loop transfer function, and Gcl is a corresponding
    closed-loop transfer function.

    Parameters
    ----------
    mags : array-like
        Array of magnitudes in dB of the M-circles
    phase_min : degrees
        Minimum phase in degrees of the N-circles
    phase_max : degrees
        Maximum phase in degrees of the N-circles

    Returns
    -------
    contours : complex array
        Array of complex numbers corresponding to the contours.
    """
    # Convert magnitudes and phase range into a grid suitable for
    # building contours
    phases = sp.radians(sp.linspace(phase_min, phase_max, 2000))
    Gcl_mags, Gcl_phases = sp.meshgrid(10.0**(mags/20.0), phases)
    return closed_loop_contours(Gcl_mags, Gcl_phases) 

def star_graph_lumped(N, tau, gamma, I0, tmin, tmax, tcount):
    times = scipy.linspace(tmin, tmax, tcount)
    #    [[central node infected] + [central node susceptible]]
    #X = [Y_1^1, Y_1^2, ..., Y_1^{N}, Y_2^0, Y_2^1, ..., Y_2^{N-1}]
    X0 = scipy.zeros(2*N)  #length 2*N of just 0 entries
    X0[I0]=I0*1./N #central infected, + I0-1 periph infected prob
    X0[N+I0] = 1-I0*1./N #central suscept + I0 periph infected
    X = EoN.my_odeint(star_graph_dX, X0, times, args = (tau, gamma, N))
    #X looks like [[central susceptible,k periph] [ central inf, k-1 periph]] x T

    central_inf = X[:,:N]
    central_susc = X[:,N:]

    I = scipy.array([ sum(k*central_susc[t][k] for k in range(N))
            + sum((k+1)*central_inf[t][k] for k in range(N))
            for t in range(len(X))])
    S = N-I
    return times, S, I 

def process_degree_distribution(N, Pk, color, Psi, DPsi, symbol, label, count):
    report_times = scipy.linspace(0,30,3000)
    sums = 0*report_times
    for cnt in range(count):
        G = generate_network(Pk, N)
        t, S, I, R = EoN.fast_SIR(G, tau, gamma, rho=rho)
        plt.plot(t, I*1./N, '-', color = color, 
                                alpha = 0.1, linewidth=1)
        subsampled_I = EoN.subsample(report_times, t, I)
        sums += subsampled_I*1./N
    ave = sums/count
    plt.plot(report_times, ave, color = 'k')
    
    #Do EBCM    
    N= G.order()#N is arbitrary, but included because our implementation of EBCM assumes N is given.
    t, S, I, R = EoN.EBCM_uniform_introduction(N, Psi, DPsi, tau, gamma, rho, tmin=0, tmax=10, tcount = 41)
    plt.plot(t, I/N, symbol, color = color, markeredgecolor='k', label=label)

    for cnt in range(3):  #do 3 highlighted simulations
        G = generate_network(Pk, N)
        t, S, I, R = EoN.fast_SIR(G, tau, gamma, rho=rho)
        plt.plot(t, I*1./N, '-', color = 'k', linewidth=0.1) 

def star_graph_lumped(N, tau, gamma, I0, tmin, tmax, tcount):
    times = scipy.linspace(tmin, tmax, tcount)
    #    [[central node infected] + [central node susceptible]]
    #X = [Y_1^1, Y_1^2, ..., Y_1^{N}, Y_2^0, Y_2^1, ..., Y_2^{N-1}]
    X0 = scipy.zeros(2*N)  #length 2*N of just 0 entries
    #X0[I0]=I0*1./N #central infected, + I0-1 periph infected prob
    X0[N+I0] = 1#-I0*1./N #central suscept + I0 periph infected
    X = EoN.my_odeint(star_graph_dX, X0, times, args = (tau, gamma, N))
    #X looks like [[central susceptible,k periph] [ central inf, k-1 periph]] x T

    central_susc = X[:,N:]
    central_inf = X[:,:N]
    print(central_susc[-1][:])
    print(central_inf[-1][:])
    I = scipy.array([ sum(k*central_susc[t][k] for k in range(N))
              + sum((k+1)*central_inf[t][k] for k in range(N))
              for t in range(len(X))])
    S = N-I
    return times, S, I 

def generate_data(N, S, L):

        # generate genetics
        G = 1.0 * (sp.rand(N, S) < 0.2)
        G -= G.mean(0)
        G /= G.std(0) * sp.sqrt(G.shape[1])

        # generate latent phenotypes
        Zg = sp.dot(G, sp.randn(G.shape[1], L))
        Zn = sp.randn(N, L)

        # generate variance exapleind
        vg = sp.linspace(0.8, 0, L)

        # rescale and sum
        Zg *= sp.sqrt(vg / Zg.var(0))
        Zn *= sp.sqrt((1 - vg) / Zn.var(0))
        Z = Zg + Zn

        return Z, G 

def _draw_vertical_guides_(self, canvas, axis_color=pyx.color.cmyk.Gray):
        """
        draws vertical guides
        """
        p = self.grid_box.params
        if p['vertical_guides']:
            line = pyx.path.path()
            nr = p['vertical_guide_nr']
            for x in scipy.linspace(self.grid_box.x_left, self.grid_box.x_right, nr):
                xt1 = self._give_trafo_x_(x, self.grid_box.y_top)
                yt1 = self._give_trafo_y_(x, self.grid_box.y_top)
                xt2 = self._give_trafo_x_(x, self.grid_box.y_bottom)
                yt2 = self._give_trafo_y_(x, self.grid_box.y_bottom)
                line.append(pyx.path.moveto(xt1, yt1))
                line.append(pyx.path.lineto(xt2, yt2))
            canvas.stroke(line, [pyx.style.linewidth.normal, pyx.style.linestyle.dotted,
                                 p['wd_axis_color']])
            # take handle
            self.ref_block_lines.append(line) 

def _draw_horizontal_guides_(self, canvas, axis_color=pyx.color.cmyk.Gray):
        """
        draws horizontal guides
        """
        p = self.grid_box.params
        if p['horizontal_guides']:
            line = pyx.path.path()
            nr = p['horizontal_guide_nr']
            for y in scipy.linspace(self.grid_box.y_top, self.grid_box.y_bottom, nr):
                xt1 = self._give_trafo_x_(self.grid_box.x_left, y)
                yt1 = self._give_trafo_y_(self.grid_box.x_left, y)
                xt2 = self._give_trafo_x_(self.grid_box.x_right, y)
                yt2 = self._give_trafo_y_(self.grid_box.x_right, y)
                line.append(pyx.path.moveto(xt1, yt1))
                line.append(pyx.path.lineto(xt2, yt2))
            canvas.stroke(line, [pyx.style.linewidth.normal, pyx.style.linestyle.dotted,
                                 p['u_axis_color']])
            self.ref_block_lines.append(line) 

def makehistogram(indata, histlen, binsize=None, therange=None):
    """

    Parameters
    ----------
    indata
    histlen
    binsize
    therange

    Returns
    -------

    """
    if therange is None:
        therange = [indata.min(), indata.max()]
    if histlen is None and binsize is None:
        thebins = 10
    elif binsize is not None:
        thebins = sp.linspace(therange[0], therange[1], (therange[1] - therange[0]) / binsize + 1, endpoint=True)
    else:
        thebins = histlen
    thehist = np.histogram(indata, thebins, therange)
    return thehist 

def upsample(inputdata, Fs_init, Fs_higher, method='univariate', intfac=False, debug=False):
    starttime = time.time()
    if Fs_higher <= Fs_init:
        print('upsample: target frequency must be higher than initial frequency')
        sys.exit()

    # upsample
    orig_x = sp.linspace(0.0, (1.0 / Fs_init) * len(inputdata), num=len(inputdata), endpoint=False)
    endpoint = orig_x[-1] - orig_x[0]
    ts_higher = 1.0 / Fs_higher
    numresamppts = int(endpoint // ts_higher + 1)
    if intfac:
        numresamppts = int(Fs_higher // Fs_init) * len(inputdata)
    else:
        numresamppts = int(endpoint // ts_higher + 1)
    upsampled_x = np.arange(0.0, ts_higher * numresamppts, ts_higher)
    upsampled_y = doresample(orig_x, inputdata, upsampled_x, method=method)
    initfilter = tide_filt.noncausalfilter(filtertype='arb', usebutterworth=False, debug=debug)
    stopfreq = np.min([1.1 * Fs_init / 2.0, Fs_higher / 2.0])
    initfilter.setfreqs(0.0, 0.0, Fs_init / 2.0, stopfreq)
    upsampled_y = initfilter.apply(Fs_higher, upsampled_y)
    if debug:
        print('upsampling took', time.time() - starttime, 'seconds')
    return upsampled_y 

def test_estimate_SIR_prob_size(self):
        print('testing estimate_SIR_prob_size')
        N = 1000000
        G = nx.fast_gnp_random_graph(N, 5. / N)
        for p in scipy.linspace(0.1, 0.5, 5):
            P, A = EoN.estimate_SIR_prob_size(G, p)
            gamma = 1.
            tau = p * gamma / (1 - p)
            P2, A2 = EoN.estimate_directed_SIR_prob_size(G, tau, 1.0)
            t, S, I, R = EoN.EBCM_discrete_from_graph(G, p)
            print("should all be approximately the same: ", R[-1] / G.order(), A, A2) 

def complete_graph_lumped(N, tau, gamma, I0, tmin, tmax, tcount):
    times = scipy.linspace(tmin, tmax, tcount)
    X0 = scipy.zeros(N+1)  #length N+1 of just 0 entries
    X0[I0]=1. #start with 100 infected.
    X = integrate.odeint(complete_graph_dX, X0, times, args = (tau, gamma, N))
    #X[t] is array whose kth entry is p(k infected| time=t).
    I = scipy.array([sum(k*Pkt[k] for k in range(len(Pkt))) for Pkt in X])
    S = N-I
    return times, S, I 

def complete_graph_lumped(N, I0, tmin, tmax, tcount):
    times = scipy.linspace(tmin, tmax, tcount)
    X0 = scipy.zeros(N+1)  #length N+1 of just 0 entries
    X0[I0]=1. #start with 100 infected.
    X = integrate.odeint(complete_graph_dX, X0, times, args = (tau, gamma, N))
    #X[t] is array whose kth entry is p(k infected| time=t).
    I = scipy.array([sum(k*Pkt[k] for k in range(len(Pkt))) for Pkt in X])
    S = N-I
    return times, S, I 

def simulate_process(graph_function, iterations, tmax, tcount, rho, kave, tau, gamma, symbol):
    Isum = scipy.zeros(tcount)
    report_times = scipy.linspace(0,tmax,tcount)
    for counter in range(iterations):
        G = graph_function()
        t, S, I = EoN.fast_SIS(G, tau, gamma, rho=rho, tmax=tmax)
        I = EoN.subsample(report_times, t, I)
        Isum += I
    plt.plot(report_times, Isum*1./(N*iterations), symbol)

#regular 

def complete_graph_lumped(N, tau, gamma, I0, tmin, tmax, tcount):
    times = scipy.linspace(tmin, tmax, tcount)
    X0 = scipy.zeros(N+1)  #length N+1 of just 0 entries
    X0[I0]=1. #start with 100 infected.
    X = integrate.odeint(complete_graph_dX, X0, times, args = (tau, gamma, N))
    #X[t] is array whose kth entry is p(k infected| time=t).
    I = scipy.array([sum(k*Pkt[k] for k in range(len(Pkt))) for Pkt in X])
    S = N-I
    return times, S, I 

def makeSpectre(transitions, sigma, step):
    """ Build a spectrum from transitions energies. For each transitions a gaussian
    function of width sigma is added in order to mimick natural broadening.
    
        :param transitions: list of transitions for readTransitions()
        :type transititions: list
        :param sigma: gaussian width in eV
        :type sigma: float
        :param step: number of absissa value
        :type step: int
        :return: absissa and spectrum value in this order
        :rtype: list, list
    
    """

    # max and min transition energies
    minval = min([val[0] for val in transitions]) - 5.0 * sigma
    maxval = max([val[0] for val in transitions]) + 5.0 * sigma
 
    # points
    npts   = int((maxval - minval) / step) + 1

    # absice
    eneval = sp.linspace(minval, maxval, npts)

    spectre = sp.zeros(npts)
    for trans in transitions:
        spectre += trans[2] * normpdf(eneval, trans[0], sigma)

    return eneval, spectre 

def n_circles(phases, mag_min=-40.0, mag_max=12.0):
    """Constant-phase contours of the function Gcl = Gol/(1+Gol), where
    Gol is an open-loop transfer function, and Gcl is a corresponding
    closed-loop transfer function.

    Parameters
    ----------
    phases : array-like
        Array of phases in degrees of the N-circles
    mag_min : dB
        Minimum magnitude in dB of the N-circles
    mag_max : dB
        Maximum magnitude in dB of the N-circles

    Returns
    -------
    contours : complex array
        Array of complex numbers corresponding to the contours.
    """
    # Convert phases and magnitude range into a grid suitable for
    # building contours
    mags = sp.linspace(10**(mag_min/20.0), 10**(mag_max/20.0), 2000)
    Gcl_phases, Gcl_mags = sp.meshgrid(sp.radians(phases), mags)
    return closed_loop_contours(Gcl_mags, Gcl_phases)


# Function aliases 

def plot(self, r=None):
        """Plot the current DIPPR correlation equation
        r correlation coefficient, optional parameter for reuse the code in the
        fit data procedure
        """
        array = self.value
        t = list(linspace(array[-2], array[-1], 100))
        kw = {}
        kw["Tc"] = self.parent.Tc.value
        kw["M"] = self.parent.M.value
        var = [DIPPR(self.prop, ti, array[:-2], **kw) for ti in t]
        dialog = Plot()
        dialog.addData(t, var)
        if self.data:
            dialog.addData(self.t, self.data, "ro")
        dialog.plot.ax.grid(True)
        dialog.plot.ax.set_title(self.title(), size="14")

        # Annotate the equation
        formula = self.formula_DIPPR(array[0], array[1:-2])
        dialog.plot.ax.annotate(formula, (0.05, 0.9), xycoords='axes fraction',
                                size="10", va="center")

        # Add correlation coefficient
        if r:
            dialog.plot.ax.annotate(
                "$r^2=%0.5f$" % r, (0.05, 0.8), xycoords='axes fraction',
                size="10", va="center")

        dialog.plot.ax.set_xlabel("T, K", ha='right', size="12")
        ylabel = "$"+self.proptex+r",\;"+self.unit.text()+"$"
        dialog.plot.ax.set_ylabel(ylabel, ha='right', size="12")
        dialog.exec_() 

def plotSoundWave(rate, sample):
    """
        Plots a given sound wave.
    """

    t = scipy.linspace(0, 2, 2 * rate, endpoint=False)
    pylab.figure('Sound wave')
    T = int(0.0001 * rate)
    pylab.plot(t[:T], sample[:T],)
    pylab.show() 

def plotSoundWave(rate, sample):
    """
        Plots a given sound wave.
    """

    t = scipy.linspace(0, 2, 2*rate, endpoint=False)
    pylab.figure('Sound wave')
    T = int(0.0001*rate)
    pylab.plot(t[:T], sample[:T],)
    pylab.show() 

def plot_models(x, y, models, fname, mx=None, ymax=None, xmin=None):

    plt.figure(num=None, figsize=(8, 6))
    plt.clf()
    plt.scatter(x, y, s=10)
    plt.title("Web traffic over the last month")
    plt.xlabel("Time")
    plt.ylabel("Hits/hour")
    plt.xticks(
        [w * 7 * 24 for w in range(10)], ['week %i' % w for w in range(10)])

    if models:
        if mx is None:
            mx = sp.linspace(0, x[-1], 1000)
        for model, style, color in zip(models, linestyles, colors):
            # print "Model:",model
            # print "Coeffs:",model.coeffs
            plt.plot(mx, model(mx), linestyle=style, linewidth=2, c=color)

        plt.legend(["d=%i" % m.order for m in models], loc="upper left")

    plt.autoscale(tight=True)
    plt.ylim(ymin=0)
    if ymax:
        plt.ylim(ymax=ymax)
    if xmin:
        plt.xlim(xmin=xmin)
    plt.grid(True, linestyle='-', color='0.75')
    plt.savefig(fname)

# first look at the data 

def plot_nullcline(self, nullc, style, lw=1, N=100, marker='', figname=None):
        if not self.do_display:
            return
        self.setup(figname)
        x_data = nullc.array[:,0]
        y_data = nullc.array[:,1]
        xs = np.sort( np.concatenate( (np.linspace(x_data[0], x_data[-1], N), x_data) ) )
        ys = nullc.spline(xs)
        pp.plot(xs, ys, style, linewidth=lw)
        if marker != '':
            pp.plot(x_data, y_data, style[0]+marker)
        self.teardown() 

def test_Gillespie_complex_contagion(self):
        def transition_rate(G, node, status, parameters):
            # this function needs to return the rate at which ``node`` changes status
            #
            r = parameters[0]
            if status[node] == 'S' and len([nbr for nbr in G.neighbors(node) if status[nbr] == 'I']) > 1:
                return 1
            else:  # status[node] might be 0 or length might be 0 or 1.
                return 0

        def transition_choice(G, node, status, parameters):
            # this function needs to return the new status of node.  We assume going
            # in that we have already calculated it is changing status.
            #
            # this function could be more elaborate if there were different
            # possible transitions that could happen.  However, for this model,
            # the 'I' nodes aren't changing status, and the 'S' ones are changing to 'I'
            # So if we're in this function, the node must be 'S' and becoming 'I'
            #
            return 'I'

        def get_influence_set(G, node, status, parameters):
            # this function needs to return any node whose rates might change
            # because ``node`` has just changed status.  That is, which nodes
            # might ``node`` influence?
            #
            # For our models the only nodes a node might affect are the susceptible neighbors.

            return {nbr for nbr in G.neighbors(node) if status[nbr] == 'S'}

        parameters = (2,)  # this is the threshold.  Note the comma.  It is needed
        # for python to realize this is a 1-tuple, not just a number.
        # ``parameters`` is sent as a tuple so we need the comma.

        N = 60000
        deg_dist = [2, 4, 6] * int(N / 3)
        G = nx.configuration_model(deg_dist)

        for rho in np.linspace(3. / 80, 7. / 80, 8):  # 8 values from 3/80 to 7/80.
            print(rho)
            IC = defaultdict(lambda: 'S')
            for node in G.nodes():
                if np.random.random() < rho:  # there are faster ways to do this random selection
                    IC[node] = 'I'

            t, S, I = EoN.Gillespie_complex_contagion(G, transition_rate, transition_choice,
                                                      get_influence_set, IC, return_statuses=('S', 'I'),
                                                      parameters=parameters)

            plt.plot(t, I)

        plt.savefig('test_Gillespie_complex_contagion') 

def  extract_hist(cls, rawfilename, winsize, nbins):
    # This function extracts the histogram features from the image

    im  = Image.open(rawfilename)
    
    (width, height) = im.size
    npixels = width * height
    pix = scipy.array(im)

    # Generate one feature vector (histogram) per pixel
    #winsize = 20  # for test.pgm
    #winsize = 0  # for RGB
    halfwin = int(winsize/2)

    bins    = scipy.linspace(0, 255, nbins)

    # Only use windows that are fully populated
    mywidth  = width-winsize
    myheight = height-winsize
    #data     = scipy.zeros((nbins-1, mywidth * myheight))
    #data     = scipy.zeros((3*winsize*winsize, mywidth * myheight))
    data    = []
    labels  = []

    # Pick up all windows, stepping by half of the window size
    for y in range(halfwin, height-halfwin, int(halfwin/2)):
      for x in range(halfwin, width-halfwin, int(halfwin/2)):
        # Read in data in row-major order
        ind = (y-halfwin)*mywidth + (x-halfwin)
        #data[:,ind] = \
        #    scipy.histogram(pix[y-halfwin:y+halfwin,
        #                        x-halfwin:x+halfwin],
        #                        bins)[0]
        # Just RGB
        #data[:,ind] = pix[y,x]
        # RGB window
        #data[:,ind] = pix[y-halfwin:y+halfwin,x-halfwin:x+halfwin].flat
        hist_features = TCData.extract_hist_subimg(pix[y-halfwin:y+halfwin,x-halfwin:x+halfwin])
        if data == []:
          data = hist_features.reshape(-1,1)
        else:
          data = scipy.concatenate((data, hist_features.reshape(-1,1)),1)
        labels    += ['(%d,%d)' % (y,x)]

    return (data, labels, width, height) 

def dotwostepresample(orig_x, orig_y, intermed_freq, final_freq, method='univariate', antialias=True, debug=False):
    """

    Parameters
    ----------
    orig_x
    orig_y
    intermed_freq
    final_freq
    method
    debug

    Returns
    -------
    resampled_y

    """
    if intermed_freq <= final_freq:
        print('intermediate frequency must be higher than final frequency')
        sys.exit()

    # upsample
    starttime = time.time()
    endpoint = orig_x[-1] - orig_x[0]
    init_freq = len(orig_x) / endpoint
    intermed_ts = 1.0 / intermed_freq
    numresamppts = int(endpoint // intermed_ts + 1)
    intermed_x = intermed_ts * sp.linspace(0.0,  1.0 * numresamppts, numresamppts, endpoint=False)
    intermed_y = doresample(orig_x, orig_y, intermed_x, method=method)
    if debug:
        print('init_freq, intermed_freq, final_freq:', init_freq, intermed_freq, final_freq)
        print('intermed_ts, numresamppts:', intermed_ts, numresamppts)
        print('upsampling took', time.time() - starttime, 'seconds')

    # antialias and ringstop filter
    if antialias:
        starttime = time.time()
        aafilterfreq = np.min([final_freq, init_freq]) / 2.0
        aafilter = tide_filt.noncausalfilter(filtertype='arb', usebutterworth=False, debug=debug)
        aafilter.setfreqs(0.0, 0.0, 0.95 * aafilterfreq, aafilterfreq)
        antialias_y = aafilter.apply(intermed_freq, intermed_y)
        if debug:
            print('antialiasing took', time.time() - starttime, 'seconds')
    else:
        antialias_y = intermed_y

    # downsample
    starttime = time.time()
    final_ts = 1.0 / final_freq
    numresamppts = np.ceil(endpoint / final_ts)
    #final_x = np.arange(0.0, final_ts * numresamppts, final_ts)
    final_x = final_ts * sp.linspace(0.0, 1.0 * numresamppts, numresamppts, endpoint=False)
    resampled_y = doresample(intermed_x, antialias_y, final_x, method=method)
    if debug:
        print('downsampling took', time.time() - starttime, 'seconds')
    return resampled_y 

def morlet(M, w=5.0, s=1.0, complete=True):
    """
    Complex Morlet wavelet.

    Parameters
    ----------
    M : int
        Length of the wavelet.
    w : float, optional
        Omega0. Default is 5
    s : float, optional
        Scaling factor, windowed from ``-s*2*pi`` to ``+s*2*pi``. Default is 1.
    complete : bool, optional
        Whether to use the complete or the standard version.

    Returns
    -------
    morlet : (M,) ndarray

    See Also
    --------
    scipy.signal.gausspulse

    Notes
    -----
    The standard version::

        pi**-0.25 * exp(1j*w*x) * exp(-0.5*(x**2))

    This commonly used wavelet is often referred to simply as the
    Morlet wavelet.  Note that this simplified version can cause
    admissibility problems at low values of `w`.

    The complete version::

        pi**-0.25 * (exp(1j*w*x) - exp(-0.5*(w**2))) * exp(-0.5*(x**2))

    This version has a correction
    term to improve admissibility. For `w` greater than 5, the
    correction term is negligible.

    Note that the energy of the return wavelet is not normalised
    according to `s`.

    The fundamental frequency of this wavelet in Hz is given
    by ``f = 2*s*w*r / M`` where `r` is the sampling rate.
    
    Note: This function was created before `cwt` and is not compatible
    with it.

    """
    x = linspace(-s * 2 * pi, s * 2 * pi, M)
    output = exp(1j * w * x)

    if complete:
        output -= exp(-0.5 * (w**2))

    output *= exp(-0.5 * (x**2)) * pi**(-0.25)

    return output 

def cwt(data, wavelet, widths):
    """
    Continuous wavelet transform.

    Performs a continuous wavelet transform on `data`,
    using the `wavelet` function. A CWT performs a convolution
    with `data` using the `wavelet` function, which is characterized
    by a width parameter and length parameter.

    Parameters
    ----------
    data : (N,) ndarray
        data on which to perform the transform.
    wavelet : function
        Wavelet function, which should take 2 arguments.
        The first argument is the number of points that the returned vector
        will have (len(wavelet(length,width)) == length).
        The second is a width parameter, defining the size of the wavelet
        (e.g. standard deviation of a gaussian). See `ricker`, which
        satisfies these requirements.
    widths : (M,) sequence
        Widths to use for transform.

    Returns
    -------
    cwt: (M, N) ndarray
        Will have shape of (len(widths), len(data)).

    Notes
    -----
    ::

        length = min(10 * width[ii], len(data))
        cwt[ii,:] = signal.convolve(data, wavelet(length,
                                    width[ii]), mode='same')

    Examples
    --------
    >>> from scipy import signal
    >>> import matplotlib.pyplot as plt
    >>> t = np.linspace(-1, 1, 200, endpoint=False)
    >>> sig  = np.cos(2 * np.pi * 7 * t) + signal.gausspulse(t - 0.4, fc=2)
    >>> widths = np.arange(1, 31)
    >>> cwtmatr = signal.cwt(sig, signal.ricker, widths)
    >>> plt.imshow(cwtmatr, extent=[-1, 1, 31, 1], cmap='PRGn', aspect='auto',
    ...            vmax=abs(cwtmatr).max(), vmin=-abs(cwtmatr).max())
    >>> plt.show()

    """
    output = np.zeros([len(widths), len(data)])
    for ind, width in enumerate(widths):
        wavelet_data = wavelet(min(10 * width, len(data)), width)
        output[ind, :] = convolve(data, wavelet_data,
                                  mode='same')
    return output 

def cwt(data, wavelet, widths):
    """
    Continuous wavelet transform.

    Performs a continuous wavelet transform on `data`,
    using the `wavelet` function. A CWT performs a convolution
    with `data` using the `wavelet` function, which is characterized
    by a width parameter and length parameter.

    Parameters
    ----------
    data : (N,) ndarray
        data on which to perform the transform.
    wavelet : function
        Wavelet function, which should take 2 arguments.
        The first argument is the number of points that the returned vector
        will have (len(wavelet(length,width)) == length).
        The second is a width parameter, defining the size of the wavelet
        (e.g. standard deviation of a gaussian). See `ricker`, which
        satisfies these requirements.
    widths : (M,) sequence
        Widths to use for transform.

    Returns
    -------
    cwt: (M, N) ndarray
        Will have shape of (len(widths), len(data)).

    Notes
    -----
    ::

        length = min(10 * width[ii], len(data))
        cwt[ii,:] = signal.convolve(data, wavelet(length,
                                    width[ii]), mode='same')

    Examples
    --------
    >>> from scipy import signal
    >>> import matplotlib.pyplot as plt
    >>> t = np.linspace(-1, 1, 200, endpoint=False)
    >>> sig  = np.cos(2 * np.pi * 7 * t) + signal.gausspulse(t - 0.4, fc=2)
    >>> widths = np.arange(1, 31)
    >>> cwtmatr = signal.cwt(sig, signal.ricker, widths)
    >>> plt.imshow(cwtmatr, extent=[-1, 1, 31, 1], cmap='PRGn', aspect='auto',
    ...            vmax=abs(cwtmatr).max(), vmin=-abs(cwtmatr).max())
    >>> plt.show()

    """
    output = np.zeros([len(widths), len(data)])
    for ind, width in enumerate(widths):
        wavelet_data = wavelet(min(10 * width, len(data)), width)
        output[ind, :] = convolve(data, wavelet_data,
                                  mode='same')
    return output 

def morlet(M, w=5.0, s=1.0, complete=True):
    """
    Complex Morlet wavelet.

    Parameters
    ----------
    M : int
        Length of the wavelet.
    w : float, optional
        Omega0. Default is 5
    s : float, optional
        Scaling factor, windowed from ``-s*2*pi`` to ``+s*2*pi``. Default is 1.
    complete : bool, optional
        Whether to use the complete or the standard version.

    Returns
    -------
    morlet : (M,) ndarray

    See Also
    --------
    scipy.signal.gausspulse

    Notes
    -----
    The standard version::

        pi**-0.25 * exp(1j*w*x) * exp(-0.5*(x**2))

    This commonly used wavelet is often referred to simply as the
    Morlet wavelet.  Note that this simplified version can cause
    admissibility problems at low values of `w`.

    The complete version::

        pi**-0.25 * (exp(1j*w*x) - exp(-0.5*(w**2))) * exp(-0.5*(x**2))

    This version has a correction
    term to improve admissibility. For `w` greater than 5, the
    correction term is negligible.

    Note that the energy of the return wavelet is not normalised
    according to `s`.

    The fundamental frequency of this wavelet in Hz is given
    by ``f = 2*s*w*r / M`` where `r` is the sampling rate.
    
    Note: This function was created before `cwt` and is not compatible
    with it.

    """
    x = linspace(-s * 2 * pi, s * 2 * pi, M)
    output = exp(1j * w * x)

    if complete:
        output -= exp(-0.5 * (w**2))

    output *= exp(-0.5 * (x**2)) * pi**(-0.25)

    return output