def filtered_cross_corr(signal1,signal2,bins,smoothing=10,timestep=1):
    """Get the cross-correlation between the signals, first filtering
    the fourier transforms (smoothed top-hat), chopping the fourier
    signals into "bins" bins"""
    signal1 -= signal1.mean()
    signal2 -= signal2.mean()
    x1 = rfft(signal1)
    x2 = rfft(signal2)
    assert len(x1)==len(x2)
    startfreq = arange(1,len(x1),len(x1)/(bins+1))
    position = arange(len(x1))
    freq = fftfreq(len(x1),timestep)
    freqout = zeros(bins)*1.
    out = zeros((bins,len(signal1)))*1.0
    att = ones(len(x1))*1.
    for i in range(bins):
        att[:startfreq[i]] = 0
        att[startfreq[i]:startfreq[i+1]] = 1
        att[startfreq[i+1]:] = 0
        freqout[i] = mean(freq*att[:len(freq)])
        att = smooth(att,smoothing)
        att[0] = 0
        x1dash = x1*att
        sig1dash = irfft(x1dash,len(signal1))
        x2dash = x2*att
        sig2dash = irfft(x2dash,len(signal2))
        out[i] = correlate(sig1dash,sig2dash,'same')
    lag = arange(-len(x2),len(x2),1.)*timestep
    return lag,freqout,out        

def filter_xyz(data): 
	xdata=[a[0] for a in data]
	ydata=[a[1] for a in data]
	zdata=[a[2] for a in data]
	samples=len(ydata)
	xdata=array(xdata)
	ydata=array(ydata)
	zdata=array(zdata)
	try:
#----------------------------------
		xfft=(fft.rfft(xdata))
		yfft=(fft.rfft(ydata))
		zfft=(fft.rfft(zdata))
	

#-------------filtering part --------------
		cutoff=samples/3*2
		xfft=ffilter(xfft,cutoff)	
		yfft=ffilter(yfft,cutoff)
		zfft=ffilter(zfft,cutoff)

		nxdata=fft.irfft(xfft)	
		nydata=fft.irfft(yfft)
		nzdata=fft.irfft(zfft)
	except:
		raise ValueError('null value')
	size=len(nxdata)
	data=[[nxdata[i],nydata[i],nzdata[i]] for i in range(size)]
	return data

 def macro_mismatch(self,p1,p2):
     """
     Performs double convolution with two different periods to calculate
     macroscopic average of a charge density along the z-axis.
     """
     from numpy.fft import rfft,irfft
     # Convert periods to direct units, if given in cartesian
     if p1 > 1.0:
         p1 = p1/self.unitcell.cell_vec[2,2]
     if p2 > 1.0:
         p2 = p2/self.unitcell.cell_vec[2,2]
     # Create xy-plane averaged density
     micro_z = self.integrate_z_density()
     # Create convolutions
     z_pos = np.linspace(0,1,len(micro_z))
     # Find index of lowest lower bound for p1
     low = 1.-p1/2.
     i1 = len(z_pos)-1
     while True:
         if z_pos[i1] <= low:
             i1 += 1
             break
         i1 -= 1
     #Find index of lowest upper bound for p1
     high = p1/2.
     j1 = 0
     while True:
         if z_pos[j1] >= high:
             j1 -= 1
             break
         j1 += 1
     # Find index of lowest lower bound for p2
     low = 1.-p2/2.
     i2 = len(z_pos)-1
     while True:
         if z_pos[i2] <= low:
             i2 += 1
             break
         i2 -= 1
     #Find index of lowest upper bound for p2
     high = p2/2.
     j2 = 0
     while True:
         if z_pos[j2] >= high:
             j2 -= 1
             break
         j2 += 1
     conv1 = np.zeros(len(micro_z))
     conv1[0:j1+1] = np.ones(j1+1)
     conv1[i1:] = np.ones(len(conv1[i1:]))
     conv2 = np.zeros(len(micro_z))
     conv2[0:j2+1] = np.ones(j2+1)
     conv2[i2:] = np.ones(len(conv2[i2:]))
     # Perform convolutions in Fourier Space
     f_micro_z = rfft(micro_z)
     f_conv1 = rfft(conv1)
     f_conv2 = rfft(conv2)
     f_macro = f_conv2*f_conv1*f_micro_z
     macro_z = irfft(f_macro)/float(np.sum(conv1))/float(np.sum(conv2))
     return macro_z

 def track_memory(self):
     '''
     Calculates the induced voltage energy kick to particles taking into
     account multi-turn induced voltage plus inductive impedance contribution.
     '''
     
     # Contribution from multi-turn induced voltage.
     self.array_memory *= np.exp(self.omegaj_array_memory * self.rev_time_array[self.counter_turn])
     induced_voltage = irfft(self.array_memory + rfft(self.slices.n_macroparticles, self.n_points_fft) * self.sum_impedances_memory, self.n_points_fft)
     self.induced_voltage = self.coefficient * induced_voltage[:self.slices.n_slices]
     induced_voltage[self.len_array_memory:]=0
     self.array_memory = rfft(induced_voltage, self.n_points_fft)
     
     # Contribution from inductive impedance
     if self.inductive_impedance_on:  
         self.induced_voltage_list[self.index_inductive_impedance].induced_voltage_generation(self.beam, 'slice_frame')
         self.induced_voltage += self.induced_voltage_list[self.index_inductive_impedance].induced_voltage
     
     # Induced voltage energy kick to particles through linear interpolation
     libfib.linear_interp_kick(self.beam.dt.ctypes.data_as(ctypes.c_void_p),
                               self.beam.dE.ctypes.data_as(ctypes.c_void_p), 
                               self.induced_voltage.ctypes.data_as(ctypes.c_void_p), 
                               self.slices.bin_centers.ctypes.data_as(ctypes.c_void_p), 
                               ctypes.c_uint(self.slices.n_slices),
                               ctypes.c_uint(self.beam.n_macroparticles),
                               ctypes.c_double(0.))
     # Counter update
     self.counter_turn += 1

def fft_lowpass(nelevation, low_bound, high_bound):
    """ Performs a low pass filter on the nelevation series.
    low_bound and high_bound specifes the boundary of the filter.
    """
    import numpy.fft as F
    if len(nelevation) % 2:
        result = F.rfft(nelevation, len(nelevation))
    else:
        result = F.rfft(nelevation)
    freq = F.fftfreq(len(nelevation))[:len(nelevation)/2]
    factor = np.ones_like(result)
    factor[freq > low_bound] = 0.0
    sl = np.logical_and(high_bound < freq, freq < low_bound)
    
    a = factor[sl]
    # Create float array of required length and reverse
    a = np.arange(len(a) + 2).astype(float)[::-1]
    # Ramp from 1 to 0 exclusive
    a = (a/a[0])[1:-1]

    # Insert ramp into factor
    factor[sl] = a
    print factor[sl]
    result = result * factor
    print result
    print 'resultnog=', len(result)
    relevation = F.irfft(result, len(nelevation))
    print 'result=', len(relevation)
    return relevation

def record2vecs(File):
    # Read the rr data (Note that it is sampled irregularly and
    # consists of the r-times in centiseconds) and return a high
    # frequency spectrogram.
    import cinc2000
    from numpy.fft import rfft, irfft
    data = []
    for line in File:
        data.append(float(line)/100)
    File.close()
    # Now data[i] is an r-time in seconds.
    hrd = cinc2000.R_times2Dev(data)
    # hrd are heart rate deviations sampled at 2 Hz
    pad = np.zeros(((Glength+len(hrd)),),np.float64)
    pad[Gl_2:len(hrd)+Gl_2] = hrd
    # Now pad is hrd with Gl_2 zeros before and after
    N_out = len(hrd)//RDt
    result = np.zeros((N_out,Fw/2),np.float64)
    mags = []
    for k in range(N_out):
        i = int(RDt*k)
        WD = Gw*pad[i:i+Glength] # Multiply data by window fuction
        FT = rfft(WD,n=Fw)
        SP = np.conjugate(FT)*FT  # Periodogram
        temp = rfft(SP,n=Fw//2)
        SP = irfft(temp[0:int(0.1*Fw/2)],n=Fw/2)
        # Low pass filter in freq domain.  Pass below 0.1 Nyquist
        temp = SP.real[:Fw/2]
        mag = math.sqrt(np.dot(temp,temp))
        result[k,:] = temp/mag
        mags.append(math.log(mag))
        # result[k,:] is a unit vector and a smoothed periodogram
    return [result,mags]

def fftfit_err(data_prof, tmp_prof, tau=0., scale=1.):
    """
    Estimates the uncertainties in the fftfit parameters.
    """

    tmp_fft = nf.rfft(tmp_prof)
    amp_tmp = np.absolute(tmp_fft)
    phs_tmp = np.angle(tmp_fft)

    dat_fft = nf.rfft(data_prof)
    amp_dat = np.absolute(dat_fft)
    phs_dat = np.angle(dat_fft)
    k = np.linspace(0, len(amp_dat)-1, num=len(amp_dat))

    diff = diffprof(data_prof, tmp_prof)
    diff_rms_sq = np.sum(diff**2) / len(diff)

    # compute Equations A10 and A11 from Taylor (1992).
    var_b = diff_rms_sq / (2 * scale * np.sum(k[1:]**2 * amp_tmp[1:] * amp_dat[1:] * \
            np.cos(phs_tmp[1:] - phs_dat[1:] + k[1:] * (tau * 2 * np.pi))))
    sig_scale = np.sqrt(var_b)
    var_t = diff_rms_sq / (2 * np.sum(amp_dat[1:]**2))
    sig_shift = np.sqrt(var_t)

    # propagate uncertainty in scale to get uncertainty in baseline. 
    sig_offset = sig_scale * amp_tmp[0] / len(data_prof)
    
    return sig_offset, sig_shift, sig_scale

def dispersion(x,t,d,fc1,fp1,fp2,fc2,alpha=1):
	from matplotlib.pyplot import ginput, plot, close, grid
	from numpy.fft import rfft,fft
	from numpy import finfo, zeros, hstack, unwrap, double, pi, angle, arctan2, imag, real
	close('all')
	eps=finfo(double).tiny
	fc1=1e6*fc1
	fp1=1e6*fp1
	fp2=1e6*fp2
	fc2=1e6*fc2
	dt=abs(t[-1]-t[-2])
	x=x-x[0]
	plot(t,x)
	grid(True)
	tx=ginput(4)
	close()
	x1=x[(t>=tx[0][0])&(t<=tx[1][0])]
	x2=x[(t>=tx[2][0])&(t<=tx[3][0])]
	toff=t[(t>=tx[1][0])&(t<=tx[2][0])]
	N1=len(x1)
	N2=len(x2)
	N3=len(toff)
	x1=hstack((x1*tukeywin(N1,alpha),zeros(N2+N3)))
	x2=hstack((zeros(N1+N3),x2*tukeywin(N2,alpha)))
	x1=bpfilter(hstack((x1*tukeywin(N1,alpha),zeros(N2+N3))),dt,fc1,fp1,fp2,fc2)
	x2=bpfilter(hstack((zeros(N1+N3),x2*tukeywin(N2,alpha))),dt,fc1,fp1,fp2,fc2)
	H=-rfft(x2)/rfft(x1)
	f=freqs(len(x1),dt)
	phi=unwrap(angle(H));
	phi=phi[(f>=fp1)&(f<=fp2)]
	f=f[(f>=fp1)&(f<=fp2)]
	c=-(4.*d*1e-3*pi*f)/phi
	f=1e-6*f[(f>=fp1)&(f<=fp2)]
	return f,c,phi

def cross_correlate(histogram, template, n_harmonics=None, upsample=16):
    """Find the shift required to align histogram with template

    """
    n = max(len(histogram),len(template))

    h_ft = rfft(histogram)
    t_ft = rfft(template)

    if len(h_ft)<len(t_ft):
        h_ft = concatenate((h_ft,zeros(len(t_ft)-len(h_ft))))
    elif len(t_ft)<len(h_ft):
        t_ft = concatenate((t_ft,zeros(len(h_ft)-len(t_ft))))

    if n_harmonics is not None:
        h_ft[n_harmonics+1:]*=0.
        t_ft[n_harmonics+1:]*=0.
    h_ft[0] = 0
    t_ft[0] = 0

    cross_correlations = irfft(conjugate(h_ft)*t_ft,n*upsample)
    shift = argmax(cross_correlations)/float(len(cross_correlations))

    assert 0<=shift<1

    #FIXME: warn if double-peaked

    return shift

def generalized_cross_correlation(d0, d1):
    # substract the means
    # (in order to get a normalized cross-correlation at the end)
    d0 -= d0.mean()
    d1 -= d1.mean()

    # Hann window to mitigate non-periodicity effects
    window = numpy.hanning(len(d0))

    # compute the cross-correlation
    D0 = rfft(d0 * window)
    D1 = rfft(d1 * window)
    D0r = D0.conjugate()
    G = D0r * D1
    # G = (G==0.)*1e-30 + (G<>0.)*G
    # W = 1. # frequency unweighted
    # W = 1./numpy.abs(G) # "PHAT"
    absG = numpy.abs(G)
    m = max(absG)
    W = 1. / (1e-10 * m + absG)
    # D1r = D1.conjugate(); G0 = D0r*D0; G1 = D1r*D1; W = numpy.abs(G)/(G0*G1) # HB weighted
    Xcorr = irfft(W * G)
    # Xcorr_unweighted = irfft(G)
    # numpy.save("d0.npy", d0)
    # numpy.save("d1.npy", d1)
    # numpy.save("Xcorr.npy", Xcorr)

    return Xcorr

    def FFT(self):
        for o in self.object_list:
            x_list = []
            y_list = []
            # print o.chaincode
            print "----"
            print o.object_name + "_1 " + str(o.point_list[0][0]) + " " + str(o.point_list[0][1]) + " 1 1"
            for d in o.chaincode:
                print d,
            print "-1"
            print "----"
            x = y = 0
            for d in o.chaincode:
                x += DIRECTION_MATRIX[d][0]
                y -= DIRECTION_MATRIX[d][1]
                x_list.append(x)
                y_list.append(y)
            n_point = len(x_list)
            x_fft_result = fft.rfft(x_list)
            y_fft_result = fft.rfft(y_list)
            for i in range(20):
                x = x_fft_result[i]
                y = y_fft_result[i]
                print "%e %e %e %e" % ( x.real, x.imag, y.real, y.imag )

            x_list_2 = fft.irfft(x_fft_result[0:10], n_point)
            y_list_2 = fft.irfft(y_fft_result[0:10], n_point)

            print "x"
            print x_list
            print x_list_2
            print "y"
            print y_list
            print y_list_2

def gaussian_convolution(data, ijk_sdev, cyclic = False, cutoff = 5,
                         task = None):

  from numpy import array, single as floatc, multiply, swapaxes
  c = array(data, floatc)

  from numpy.fft import rfft, irfft
  for axis in range(3):           # Transform one axis at a time.
    size = c.shape[axis]
    sdev = ijk_sdev[2-axis]       # Axes i,j,k are 2,1,0.
    hw = min(size/2, int(cutoff*sdev+1)) if cutoff else size/2
    nzeros = 0 if cyclic else hw  # Zero-fill for non-cyclic convolution.
    if nzeros > 0:
      # FFT performance is much better (up to 10x faster in numpy 1.2.1)
      # than other sizes.
      nzeros = efficient_fft_size(size + nzeros) - size
    g = gaussian(sdev, size + nzeros)
    g[hw:-hw] = 0
    fg = rfft(g)                  # Fourier transform of 1-d gaussian.
    cs = swapaxes(c, axis, 2)     # Make axis 2 the FT axis.
    s0 = cs.shape[0]
    for p in range(s0):  # Transform one plane at a time.
      cp = cs[p,...]
      try:
        ft = rfft(cp, n=len(g))   # Complex128 result, size n/2+1
      except ValueError, e:
        raise MemoryError, e      # Array dimensions too large.
      multiply(ft, fg, ft)
      cp[:,:] = irfft(ft)[:,:size] # Float64 result
      if task:
        pct = 100.0 * (axis + float(p)/s0) / 3.0
        task.updateStatus('%.0f%%' % pct)

 def periodic_interp(self,data, zoomfact, window='hanning', alpha=6.0):
     """
     Return a periodic, windowed, sinc-interpolation of the data which
     is oversampled by a factor of 'zoomfact'.
     """
     zoomfact = int(zoomfact)
     if (zoomfact < 1):
         #print "zoomfact must be >= 1."
         return 0.0
     elif zoomfact==1:
         return data
     
     newN = len(data)*zoomfact
     # Space out the data
     comb = Num.zeros((zoomfact, len(data)), dtype='d')
     comb[0] += data
     comb = Num.reshape(Num.transpose(comb), (newN,))
     # Compute the offsets
     xs = Num.zeros(newN, dtype='d')
     xs[:newN/2+1] = Num.arange(newN/2+1, dtype='d')/zoomfact
     xs[-newN/2:]  = xs[::-1][newN/2-1:-1]
     # Calculate the sinc times window for the kernel
     if window.lower()=="kaiser":
         win = _window_function[window](xs, len(data)/2, alpha)
     else:
         win = _window_function[window](xs, len(data)/2)
     kernel = win * self.sinc(xs)
     
     if (0):
         print "would have plotted."
     return FFT.irfft(FFT.rfft(kernel) * FFT.rfft(comb))

def processScenario(m_surface, scenario, servo_id, beat, mode=0):
	tmp = []
	m_data = []
	#l = len(m_surface)
	s = len(scenario)
	ms_index = -1	
	
	for i in range(s):	# << the augmented gestures can only be as long as the scenario length
		p_data = []
		#m_data = []
		p_data_tmp = []
		mapped = False
		r = 0
		#l = len(scenario[i])
		l = len(scenario[i])
		print "length:", l
		"""if m_surface[i]['phrase'] and m_surface[i]['notes'] != []:		# If it's a (melodic) phrase...
							# just a reminder of the contents of phrase info: 
							# 'phrase_start', 'phrase_end', 'phrase_length', 'contour_amplitudes', 'contour_energy',
							# 'contour_time', 'notes'
			p_data = numpy.multiply(m_surface[i]['contour_amplitudes'], m_surface[i]['contour_energy'])
			p_data_ = [converts(int(p)) for p in p_data]
			print "scenario[%d] = " % (i), scenario[i]
			print "m_surface[%d] = " % (i), m_surface[i]
			p_data_tmp = fft.ifft(numpy.add(fft.rfft(scenario[i], n=l), fft.rfft(p_data_, n=l)))
			
			print p_data_tmp			
			
			if m_surface[i]['contour_energy'] != []: tmp.extend([accel(m_surface[i]['contour_energy'][0], servo_id)])
			tmp.extend(p_data_tmp)
		
		elif not m_surface[i]['phrase']:				# Otherwise it's a rest (phrase).
			tmp.extend([m_surface[i]['rest_duration']*beat])
		else: continue
		"""
	
		
		while not mapped:
			ms_index += 1
			if m_surface[ms_index]['phrase'] and m_surface[ms_index]['phrase_length'] != 0 and len(m_surface[ms_index]['contour_amplitudes']) > (l/2):
				 
				print m_surface[ms_index]
				#p_data = fft.irfft(fft.rfft(scenario[i], l)*fft.rfft(m_surface[ms_index]['contour_amplitudes'], l))
				p_data = fft.irfft(numpy.multiply(fft.rfft(scenario[i], l), .02*(fft.rfft(m_surface[ms_index]['contour_amplitudes'], l)[2])))#(numpy.multiply(m_surface[ms_index]['contour_amplitudes'], m_surface[ms_index]['contour_energy']), l)))
				mapped = True
				print "pdata:",p_data
				m_data.append(m_surface[ms_index]['contour_amplitudes'])
				break
			elif m_surface[ms_index]['phrase'] == False:
				#print m_surface
				#p_data = 'rest'
				#tmp.append(m_surface[ms_index]['rest_duration'])
				r += m_surface[ms_index]['rest_duration']
				print "rest"
		tmp.append(r)		
		tmp.append(p_data)		
	
	#print tmp
	return tmp, m_data

def my_multi_convolve(signals):
    shape = sum(len(x)-1 for x in signals) + 1
    fshape = _next_regular(shape)
    res = rfft(signals[0], fshape)
    for signal in signals[1:]:
        res *= rfft(signal, fshape)
    ret = irfft(res, fshape)
    return ret[:shape]

def _bode(idat, odat, fs):
	ift = rfft(idat, fs)
	oft =  rfft(odat, fs)
	trans = (ift * conjugate(oft)) / (ift * conjugate(ift) ) 
	freq = arange(trans.shape[0]).astype(Float64)*fs/(2*trans.shape[0])
	amp=abs(trans)
	phase=arctan2(trans.imag, trans.real)
	return (freq, amp, phase)

def fourier_analysis(data, ptree=None):
    time = data['time']
    current = data['current']
    voltage = data['voltage']

    # inspect data
    n = len(time)  # normalization factor for fft
    assert len(current) == n
    assert len(voltage) == n
    d = time[1] - time[0]  # inverse of the sampling rate
    # check sampling spacing is the same everywhere
    for i in range(n - 1):
        assert isclose(time[i + 1] - time[i], d, atol=1e-10, rtol=1e-10)

    # truncate signals
    if ptree:
        steps_per_cycle = ptree.get_int('steps_per_cycle')
        cycles = ptree.get_int('cycles')
        ignore_cycles = ptree.get_int('ignore_cycles')
        assert cycles > ignore_cycles
        assert n == cycles * steps_per_cycle
        time = time[ignore_cycles * steps_per_cycle:]
        current = current[ignore_cycles * steps_per_cycle:]
        voltage = voltage[ignore_cycles * steps_per_cycle:]
    else:
        time = time[int(n / 2):]
        current = current[int(n / 2):]
        voltage = voltage[int(n / 2):]

    n = len(time)
    assert len(current) == n
    assert len(voltage) == n

    if not _is_power_of_two(n):
        warn(
            "(cycles-ignore_cycles)*steps_per_cycles is not a "
            "power of 2 (most efficient for the fourier analysis)",
            RuntimeWarning)

    # perform the actual fourrier analaysis
    fft_current = fft.rfft(current) / n
    fft_voltage = fft.rfft(voltage) / n
    fft_frequency = fft.rfftfreq(n, d)

    # find the excited harmonics
    if ptree:
        harmonics = array(ptree.get_array_int('harmonics'))
        peak_indices = harmonics * (cycles - ignore_cycles)
    else:
        mx, mn = peakdet(absolute(fft_voltage), mean(absolute(fft_current)))
        peak_indices = int(mx[:, 0])
        mx, mn = peakdet(absolute(fft_voltage), mean(absolute(fft_current)))
        assert peak_indices == mx[:, 0]

    frequency = fft_frequency[peak_indices]
    impedance = fft_voltage[peak_indices] / fft_current[peak_indices]

    return [frequency, impedance]

def rms_rfftdf(df):
      X=rfft(df['accel-x'])
      Y=rfft(df['accel-y'])
      Z=rfft(df['accel-z'])  
      
      x_rfft= rms_rfft(X, len(df))
      y_rfft= rms_rfft(Y, len(df))
      z_rfft= rms_rfft(Z, len(df))
      return pd.DataFrame({'x': [x_rfft],'y': [y_rfft], 'z':[z_rfft]})

def H1(x,y,dt,frng,uwphase=True,nfreqs='All'):
    """
    Returns a frequency vector, transfer function between x and y, and the wrapped phase of the transfer
    function
    """    
    
    from numpy import zeros,conjugate,linspace, angle, unwrap, array, arange
    from numpy.fft import rfft    
    
    # X = rfft(x,int(1/(df*dt)))
   #  Y = rfft(y,int(1/(df*dt)))
    X=rfft(x)
    
    Xc=conjugate(X)
    
    H=[]
    phi=[]
    f=linspace(0.,1/(2*dt),len(X))
    
    
    df=f[1]-f[0]
    if1=NearestValue(f,frng[0])
    if2=NearestValue(f,frng[1])
    
    if nfreqs is 'All':
        
        find=arange(if1,if2,1)
        
    else:
        
        find=linspace(if1,if2,nfreqs).astype(int)
    

    for i in range(len(y)):
        
        Y=rfft(y[i])
        HH=(Y*Xc)/(X*Xc)
        
        pphi=angle(HH)
        
        if uwphase:
            pphi=unwrap(angle(HH))
            
            
        H.append(HH[find])
        phi.append(pphi[find])
    
    f=f[find]
    
    f=array(f)
    H=array(H)
    phi=array(phi)
    
    return f,H,phi

 def to_freq(self, data_type):
     #self._data1 and self._data2 need to be defined by the object
     delta_f = 1 / (self.delta_t * self.pulses_length)
     data1_freq = rfft(self._data[self._data1])
     data2_freq = rfft(self._data[self._data2])
     freq_pulses_length = data1_freq.shape[1]
     pulses_freq = data_type(freq_pulses_length, self.pulses_nb, delta_f)
     pulses_freq._data['Valim'] = self.valim
     pulses_freq._data[self._data1] = data1_freq
     pulses_freq._data[self._data2] = data2_freq
     return pulses_freq