Python scipy.signal.resample() Examples

def cutTabs(datas, part):
	#First we uniformize the GS
	minSize = 0
	for nDim in range(len(v.eName)):
		for s in part :
			if (minSize > len(datas['gs'+s][nDim]) or minSize == 0):
				minSize = len(datas['gs'+s][nDim])
	oneF = int(minSize/9)
	#We cut all tab to reach this size
	for nDim in range(len(v.eName)):
		for s in part :
			#Gold Standard Tab
			datas['gs'+s][nDim] = cutTab(datas['gs'+s][nDim],minSize)
			#Predictions tab
			for nMod in range(len(v.desc)):
				datas[s][nDim][nMod] = cutTab(datas[s][nDim][nMod],minSize)
	return datas
#End cutTabs

#Used to resample the tab 

def _wavdata_rs(self):
        if self.fs_rs is not None:
            # Number of points in resample
            ns_rs = np.int_(np.ceil(self.ns * self.fs_rs / self.fs))
            # Do resample
            # XXX: Tried using a Hamming window as a low pass filter, but it
            #      didn't seem to make a big difference, so it's not used
            #      here.
            data_rs = resample(self.wavdata, ns_rs)
            wavpath_rs = self.wavpath.split('.')[0] + '-resample-' + str(self.fs_rs) + 'Hz.wav'
            # Write resampled data to wav file
            # Convert data from 32-bit floating point to 16-bit PCM
            data_rs_int = np.int16(data_rs * 32768)
            wavfile.write(wavpath_rs, self.fs_rs, data_rs_int)
            # XXX: Was worried that Python might continue executing code
            #      before the file write is finished, but it seems like it's
            #      not an issue.
            return wavpath_rs, data_rs, data_rs_int, ns_rs
        else:
            return None, None, None, None 

def test_basic(self):
        # Some basic tests

        # Regression test for issue #3603.
        # window.shape must equal to sig.shape[0]
        sig = np.arange(128)
        num = 256
        win = signal.get_window(('kaiser', 8.0), 160)
        assert_raises(ValueError, signal.resample, sig, num, window=win)

        # Other degenerate conditions
        assert_raises(ValueError, signal.resample_poly, sig, 'yo', 1)
        assert_raises(ValueError, signal.resample_poly, sig, 1, 0)

        # test for issue #6505 - should not modify window.shape when axis ≠ 0
        sig2 = np.tile(np.arange(160), (2,1))
        signal.resample(sig2, num, axis=-1, window=win)
        assert_(win.shape == (160,)) 

def downsample(data, oldFS, newFS):
    """
    Resample data from oldFS to newFS using the scipy 'resample' function.

    Parameters
    ----------
    data : instance of pandas.core.DataFrame
        Data to resample.
    oldFS : float
        The sampling frequency of data.
    newFS : float
        The new sampling frequency.

    Returns:

    newData : instance of pandas.DataFrame
        The downsampled dataset.
    """

    newNumSamples = int((data.shape[0] / oldFS) * newFS)
    newData = pd.DataFrame(resample(data, newNumSamples))
    return newData 

def get_traces(self, channel_ids=None, start_frame=None, end_frame=None):
        start_frame_not_sampled = int(start_frame / self.get_sampling_frequency() *
                                      self._recording.get_sampling_frequency())
        start_frame_sampled = start_frame
        end_frame_not_sampled = int(end_frame / self.get_sampling_frequency() *
                                    self._recording.get_sampling_frequency())
        end_frame_sampled = end_frame
        traces = self._recording.get_traces(start_frame=start_frame_not_sampled,
                                            end_frame=end_frame_not_sampled,
                                            channel_ids=channel_ids)
        if np.mod(self._recording.get_sampling_frequency(), self._resample_rate) == 0:
            traces_resampled = signal.decimate(traces,
                                               q=int(self._recording.get_sampling_frequency() / self._resample_rate),
                                               axis=1)
        else:
            traces_resampled = signal.resample(traces, int(end_frame_sampled - start_frame_sampled), axis=1)
        return traces_resampled.astype(self._dtype) 

def resample(recording, resample_rate):
    '''
    Resamples the recording extractor traces. If the resampling rate is multiple of the sampling rate, the faster
    scipy decimate function is used.

    Parameters
    ----------
    recording: RecordingExtractor
        The recording extractor to be resampled
    resample_rate: int or float
        The resampling frequency

    Returns
    -------
    resampled_recording: ResampleRecording
        The resample recording extractor

    '''
    return ResampleRecording(
        recording=recording,
        resample_rate=resample_rate
    ) 

def preprocess(data, config):
    sr = config.sample_rate
    if sr == None:
      sr = 300
    data = np.nan_to_num(data) # removing NaNs and Infs
    from scipy.signal import resample
    data = resample(data, int(len(data) * 360 / sr) ) # resample to match the data sampling rate 360(mit), 300(cinc)
    from sklearn import preprocessing
    data = preprocessing.scale(data)
    from scipy.signal import find_peaks
    peaks, _ = find_peaks(data, distance=150)
    data = data.reshape(1,len(data))
    data = np.expand_dims(data, axis=2) # required by Keras
    return data, peaks

# predict 

def __call__(self, seq):
        if np.random.randint(2):
            return seq
        else:
            seq_aug = np.zeros(seq.shape)
            len = seq.shape[1]
            length = int(len * (1 + (random.random()-0.5)*self.sigma))
            for i in range(seq.shape[0]):
                y = resample(seq[i, :], length)
                if length < len:
                    if random.random() < 0.5:
                        seq_aug[i, :length] = y
                    else:
                        seq_aug[i, len-length:] = y
                else:
                    if random.random() < 0.5:
                        seq_aug[i, :] = y[:len]
                    else:
                        seq_aug[i, :] = y[length-len:]
            return seq_aug 

def downsampleNP(data, oldFS, newFS):
    """
    Resample data from oldFS to newFS using the scipy 'resample' function.

    Parameters
    ----------
    data : array-like
        Data to resample.
    oldFS : float
        The sampling frequency of data.
    newFS : float
        The new sampling frequency.

    Returns:

    newData : instance of pandas.DataFrame
        The downsampled dataset.
    """

    newNumSamples = int((data.shape[0] / oldFS) * newFS)
    newData = resample(data, newNumSamples)
    return newData 

def simple_resample(data, original_fs, target_fs):
    """Wrap scipy.signal.resample with a simpler API
    """
    return resample(data, num=int(data.shape[-1] * target_fs / original_fs), axis=-1) 

def apply(self, data):
        axis = data.ndim - 1
        if data.shape[-1] > self.f:
            return resample(data, self.f, axis=axis)
        return data 

def resamplingTab(tab, size):
	if (len(tab) != size):
		s = signal.resample(tab,size)
		return s
	else :
		return tab
#End resamplingTab

#Calculus of CCC 

def unimodalPredPrep(wSize, wStep, nMod):
	feats = {}
	#We need the number of line for a wStep of v.tsp
	trainLen = len(arff.load(open(v.descNorm[nMod]+"train_"+str(wSize)+"_"+str(v.tsp)+".arff","rb"))['data'])
	#We open corresponding files
	for s in v.part:	
		feats[s] = arff.load(open(v.descNorm[nMod]+s+"_"+str(wSize)+"_"+str(wStep)+".arff","rb"))
		#We put to 0 NaN values
		feats[s] = arffNan(feats[s])
		#We transform it in array
		feats[s] = np.array(feats[s]['data'])
		#We resample it to be at a wSize of v.tsp
		feats[s] = resamplingTab(feats[s], trainLen)
	return feats, trainLen
#End unimodalPredPrep 

def deweightmem(input_data, weight_fun=None, oversample_rate=1, dim=1):
    """DEWEIGHTMEM Make complex SAR uniformly weighted in one dimension.

          Parameter name    Description

          input_data        Array of complex values to deweight
          weight_fun        Description of weighting applied.  Either a
                               function that takes a single argument (number
                               of elements) and produces the weighting to
                               apply, or a vector that is the weighting
                               function sampled.
          oversample_rate   Amount of sampling beyond the ImpRespBW in the
                               processing dimension. (Default is Nyquist
                               sampling = 1).
          dim               Dimension over which to perform deweighting.
                               Default is 1.

    This implementation assumes that the data has already been "deskewed" and
    that the frequency support is centered.
    """
    # TODO: Test this function

    # Weighting only valid across ImpRespBW
    weight_size = round(input_data.shape[dim]/oversample_rate)
    if weight_fun is None:  # No weighting passed in.  Do nothing.
        return input_data
    elif callable(weight_fun):
        weighting = weight_fun(weight_size)
    elif np.array(weight_fun).ndim == 1:
        import scipy.signal as sig
        weighting = sig.resample(weight_fun, weight_size)

    weight_zp = np.ones(input_data.shape[dim])  # Don't scale outside of ImpRespBW
    weight_zp[np.floor((input_data.shape[dim]-weight_size)/2)+np.arange(weight_size)] = weighting

    # Divide out weighting in spatial frequency domain
    output_data = np.fft.fftshift(np.fft.fft(input_data, axis=dim), axes=dim)
    output_data = output_data/weight_zp
    output_data = np.fft.ifft(np.fft.ifftshift(output_data, axes=dim), axis=dim)

    return output_data 

def add_noise(config):
    noises = dict()
    noises["trainset"] = list()
    noises["testset"] = list() 
    import csv
    try:
        testlabel = list(csv.reader(open('training2017/REFERENCE.csv')))
    except:
        cmd = "curl -O https://archive.physionet.org/challenge/2017/training2017.zip"
        os.system(cmd)
        os.system("unzip training2017.zip")
        testlabel = list(csv.reader(open('training2017/REFERENCE.csv')))
    for i, label in enumerate(testlabel):
      if label[1] == '~':
        filename = 'training2017/'+ label[0] + '.mat'
        from scipy.io import loadmat
        noise = loadmat(filename)
        noise = noise['val']
        _, size = noise.shape
        noise = noise.reshape(size,)
        noise = np.nan_to_num(noise) # removing NaNs and Infs
        from scipy.signal import resample
        noise= resample(noise, int(len(noise) * 360 / 300) ) # resample to match the data sampling rate 360(mit), 300(cinc)
        from sklearn import preprocessing
        noise = preprocessing.scale(noise)
        noise = noise/1000*6 # rough normalize, to be improved 
        from scipy.signal import find_peaks
        peaks, _ = find_peaks(noise, distance=150)
        choices = 10 # 256*10 from 9000
        picked_peaks = np.random.choice(peaks, choices, replace=False)
        for j, peak in enumerate(picked_peaks):
          if peak > config.input_size//2 and peak < len(noise) - config.input_size//2:
              start,end  = peak-config.input_size//2, peak+config.input_size//2
              if i > len(testlabel)/6:
                noises["trainset"].append(noise[start:end].tolist())
              else:
                noises["testset"].append(noise[start:end].tolist())
    return noises 

def __init__(self, timeaxis, timecourse, padvalue=30.0, upsampleratio=100, doplot=False, debug=False,
                 method='univariate'):
        self.upsampleratio = upsampleratio
        self.padvalue = padvalue
        self.initstep = timeaxis[1] - timeaxis[0]
        self.initstart = timeaxis[0]
        self.initend = timeaxis[-1]
        self.hiresstep = self.initstep / np.float64(self.upsampleratio)
        self.hires_x = np.arange(timeaxis[0] - self.padvalue, self.initstep * len(timeaxis) + self.padvalue,
                                 self.hiresstep)
        self.hiresstart = self.hires_x[0]
        self.hiresend = self.hires_x[-1]
        if method == 'poly':
            self.hires_y = 0.0 * self.hires_x
            self.hires_y[int(self.padvalue // self.hiresstep) + 1:-(int(self.padvalue // self.hiresstep) + 1)] = \
                signal.resample_poly(timecourse, np.int(self.upsampleratio * 10), 10)
        elif method == 'fourier':
            self.hires_y = 0.0 * self.hires_x
            self.hires_y[int(self.padvalue // self.hiresstep) + 1:-(int(self.padvalue // self.hiresstep) + 1)] = \
                signal.resample(timecourse, self.upsampleratio * len(timeaxis))
        else:
            self.hires_y = doresample(timeaxis, timecourse, self.hires_x, method=method)
        self.hires_y[:int(self.padvalue // self.hiresstep)] = self.hires_y[int(self.padvalue // self.hiresstep)]
        self.hires_y[-int(self.padvalue // self.hiresstep):] = self.hires_y[-int(self.padvalue // self.hiresstep)]
        if debug:
            print('fastresampler __init__:')
            print('    padvalue:, ', self.padvalue)
            print('    initstep, hiresstep:', self.initstep, self.hiresstep)
            print('    initial axis limits:', self.initstart, self.initend)
            print('    hires axis limits:', self.hiresstart, self.hiresend)

        # self.hires_y[:int(self.padvalue // self.hiresstep)] = 0.0
        # self.hires_y[-int(self.padvalue // self.hiresstep):] = 0.0
        if doplot:
            fig = pl.figure()
            ax = fig.add_subplot(111)
            ax.set_title('fastresampler initial timecourses')
            pl.plot(timeaxis, timecourse, self.hires_x, self.hires_y)
            pl.legend(('input', 'hires'))
            pl.show() 

def apply(self, data):
        axis = data.ndim - 1
        out = resample(data, self.f, axis=axis, window=hann(M=data.shape[axis]))
        return out 

def get_name(self):
        return "resample%d" % self.f 

def apply(self, data):
        axis = data.ndim - 1
        if data.shape[-1] > self.f:
            return resample(data, self.f, axis=axis)
        return data 

def get_name(self):
        return "resample%dhanning" % self.f 

def init_sound(self):
        """
        Load the wavfile with :mod:`scipy.io.wavfile` ,
        converting int to float as needed.

        Create a sound table, resampling sound if needed.


        """

        fs, audio = wavfile.read(self.path)
        if audio.dtype in ['int16', 'int32']:
            audio = int_to_float(audio)

        # load file to sound table
        if self.server_type == 'pyo':
            self.dtable = pyo.DataTable(size=audio.shape[0], chnls=prefs.NCHANNELS, init=audio.tolist())

            # get server to determine sampling rate modification and duration
            server_fs = self.dtable.getServer().getSamplingRate()
            self.duration = float(self.dtable.getSize()) / float(fs)
            self.table = pyo.TableRead(table=self.dtable, freq=float(fs) / server_fs,
                                       loop=False, mul=self.amplitude)

        elif self.server_type == 'jack':
            # attenuate amplitude
            audio = audio*self.amplitude
            self.duration = float(audio.shape[0]) / fs
            # resample to match our audio server's sampling rate
            if fs != self.fs:
                new_samples = self.duration*self.fs
                audio = resample(audio, new_samples)

            self.table = audio

        self.initialized = True 

def get_name(self):
        return "resample%d" % self.f 

def to_envelopes(path,num_bands,freq_lims,downsample=True):
    """Generate amplitude envelopes from a full path to a .wav, following
    Lewandowski (2012).

    Parameters
    ----------
    filename : str
        Full path to .wav file to process.
    num_bands : int
        Number of frequency bands to use.
    freq_lims : tuple
        Minimum and maximum frequencies in Hertz to use.
    downsample : bool
        If true, envelopes will be downsampled to 120 Hz

    Returns
    -------
    2D array
        Amplitude envelopes over time
    """
    sr, proc = preproc(path,alpha=0.97)

    #proc = proc / 32768 #hack!! for 16-bit pcm
    proc = proc/sqrt(mean(proc**2))*0.03;
    bandLo = [ freq_lims[0]*exp(log(freq_lims[1]/freq_lims[0])/num_bands)**x for x in range(num_bands)]
    bandHi = [ freq_lims[0]*exp(log(freq_lims[1]/freq_lims[0])/num_bands)**(x+1) for x in range(num_bands)]

    envelopes = []
    for i in range(num_bands):
        b, a = butter(2,(bandLo[i]/(sr/2),bandHi[i]/(sr/2)), btype = 'bandpass')
        env = filtfilt(b,a,proc)
        env = abs(hilbert(env))
        if downsample:
            env = resample(env,int(ceil(len(env)/int(ceil(sr/120)))))
            #env = decimate(env,int(ceil(sr/120)))
        envelopes.append(env)
    return array(envelopes).T 

def preproc(path,sr=None,alpha=0.97):
    """Preprocess a .wav file for later processing.

    Parameters
    ----------
    path : str
        Full path to .wav file to load.
    sr : int, optional
        Sampling rate to resample at, if specified.
    alpha : float, optional
        Alpha for preemphasis, defaults to 0.97.

    Returns
    -------
    int
        Sampling rate.
    array
        Processed PCM.

    """
    oldsr,sig = wavfile.read(path)

    try:
        sig = sig[:,0]
    except IndexError:
        pass

    if sr is not None and sr != oldsr:
        t = len(sig)/oldsr
        numsamp = t * sr
        proc = resample(sig,numsamp)
    else:
        sr = oldsr
        proc = sig
    if alpha != 0:
        proc = lfilter([1., -alpha],1,proc)
    return sr,proc 

def transform(sig, train=False):
    # 前置不可或缺的步骤
    sig = resample(sig, config.target_point_num)
    # # 数据增强
    if train:
        if np.random.randn() > 0.5: sig = scaling(sig)
        if np.random.randn() > 0.5: sig = verflip(sig)
        if np.random.randn() > 0.5: sig = shift(sig)
    # 后置不可或缺的步骤
    sig = sig.transpose()
    sig = torch.tensor(sig.copy(), dtype=torch.float)
    return sig 

def resample(sig, target_point_num=None):
    '''
    对原始信号进行重采样
    :param sig: 原始信号
    :param target_point_num:目标型号点数
    :return: 重采样的信号
    '''
    sig = signal.resample(sig, target_point_num) if target_point_num else sig
    return sig 

def resample_dataset(x, rate):
    new_x = np.zeros(shape=(x.shape[0], rate))
    from scipy import signal
    for i in range(x.shape[0]):
        f = signal.resample(x[0], rate)
        new_x[i] = f
    return new_x 

def upsample_resample(wf, upsample_factor):
    wf = wf.T
    waveform_len, n_spikes = wf.shape
    traces = np.zeros((n_spikes, (waveform_len-1)*upsample_factor+1),'float32')
    for j in range(wf.shape[1]):
        traces[j] = signal.resample(wf[:,j],(waveform_len-1)*upsample_factor+1)
    return traces 

def test_array_as_window(self):
        # github issue 3603
        osfactor = 128
        sig = np.arange(128)

        win = windows.get_window(('kaiser', 8.0), osfactor // 2)
        assert_raises(ValueError, resample,
                      (sig, len(sig) * osfactor), {'window': win}) 

def __call__(self, pkg):
        pkg = format_package(pkg)
        wav = pkg['chunk']
        wav = wav.data.numpy().reshape(-1).astype(np.float32)
        warp_factor = random.random() * (self.factor_range[1] - \
                                         self.factor_range[0]) + \
                      self.factor_range[0]
        samp_warp = wav.shape[0] + int(warp_factor * wav.shape[0])
        rwav = signal.resample(wav, samp_warp)
        if len(rwav) > len(wav):
            mid_i = (len(rwav) // 2) - len(wav) // 2
            rwav = rwav[mid_i:mid_i + len(wav)]
        if len(rwav) < len(wav):
            diff = len(wav) - len(rwav)
            P = (len(wav) - len(rwav)) // 2
            if diff % 2 == 0:
                rwav = np.concatenate((np.zeros(P, ),
                                       wav,
                                       np.zeros(P, )),
                                      axis=0)
            else:
                rwav = np.concatenate((np.zeros(P, ),
                                       wav,
                                       np.zeros(P + 1, )),
                                      axis=0)
        if self.report:
            if 'report' not in pkg:
                pkg['report'] = {}
            pkg['report']['warp_factor'] = warp_factor
        pkg['chunk'] = torch.FloatTensor(rwav)
        return pkg