paulstretch_python

PaulStretch python version
Log | Files | Refs | README

paulstretch_newmethod.py (6828B)


      1 #!/usr/bin/env python
      2 #
      3 # Paul's Extreme Sound Stretch (Paulstretch) - Python version
      4 # using a new method
      5 #
      6 # by Nasca Octavian PAUL, Targu Mures, Romania
      7 # http://www.paulnasca.com/
      8 #
      9 # http://hypermammut.sourceforge.net/paulstretch/
     10 #
     11 # this file is released under Public Domain
     12 #
     13 
     14 
     15 import sys
     16 from numpy import *
     17 import scipy.io.wavfile
     18 import wave
     19 from optparse import OptionParser
     20 
     21 plot_onsets=False
     22 if plot_onsets:
     23     import matplotlib.pyplot as plt
     24 
     25 
     26 def load_wav(filename):
     27     try:
     28         wavedata=scipy.io.wavfile.read(filename)
     29         samplerate=int(wavedata[0])
     30         smp=wavedata[1]*(1.0/32768.0)
     31         smp=smp.transpose()
     32         if len(smp.shape)==1: #convert to stereo
     33             smp=tile(smp,(2,1))
     34         return (samplerate,smp)
     35     except:
     36         print ("Error loading wav: "+filename)
     37         return None
     38 
     39 
     40 
     41 def optimize_windowsize(n):
     42     orig_n=n
     43     while True:
     44         n=orig_n
     45         while (n%2)==0:
     46             n/=2
     47         while (n%3)==0:
     48             n/=3
     49         while (n%5)==0:
     50             n/=5
     51 
     52         if n<2:
     53             break
     54         orig_n+=1
     55     return orig_n
     56 
     57 def paulstretch(samplerate,smp,stretch,windowsize_seconds,onset_level,outfilename):
     58 
     59     if plot_onsets:
     60         onsets=[]
     61 
     62     nchannels=smp.shape[0]
     63 
     64     outfile=wave.open(outfilename,"wb")
     65     outfile.setsampwidth(2)
     66     outfile.setframerate(samplerate)
     67     outfile.setnchannels(nchannels)
     68 
     69     #make sure that windowsize is even and larger than 16
     70     windowsize=int(windowsize_seconds*samplerate)
     71     if windowsize<16:
     72         windowsize=16
     73     windowsize=optimize_windowsize(windowsize)
     74     windowsize=int(windowsize/2)*2
     75     half_windowsize=int(windowsize/2)
     76 
     77     #correct the end of the smp
     78     nsamples=smp.shape[1]
     79     end_size=int(samplerate*0.05)
     80     if end_size<16:
     81         end_size=16
     82 
     83     smp[:,nsamples-end_size:nsamples]*=linspace(1,0,end_size)
     84 
     85     
     86     #compute the displacement inside the input file
     87     start_pos=0.0
     88     displace_pos=windowsize*0.5
     89 
     90     #create Hann window
     91     window=0.5-cos(arange(windowsize,dtype='float')*2.0*pi/(windowsize-1))*0.5
     92 
     93     old_windowed_buf=zeros((2,windowsize))
     94     hinv_sqrt2=(1+sqrt(0.5))*0.5
     95     hinv_buf=2.0*(hinv_sqrt2-(1.0-hinv_sqrt2)*cos(arange(half_windowsize,dtype='float')*2.0*pi/half_windowsize))/hinv_sqrt2
     96 
     97     freqs=zeros((2,half_windowsize+1))
     98     old_freqs=freqs
     99 
    100     num_bins_scaled_freq=32
    101     freqs_scaled=zeros(num_bins_scaled_freq)
    102     old_freqs_scaled=freqs_scaled
    103 
    104     displace_tick=0.0
    105     displace_tick_increase=1.0/stretch
    106     if displace_tick_increase>1.0:
    107         displace_tick_increase=1.0
    108     extra_onset_time_credit=0.0
    109     get_next_buf=True
    110     while True:
    111         if get_next_buf:
    112             old_freqs=freqs
    113             old_freqs_scaled=freqs_scaled
    114 
    115             #get the windowed buffer
    116             istart_pos=int(floor(start_pos))
    117             buf=smp[:,istart_pos:istart_pos+windowsize]
    118             if buf.shape[1]<windowsize:
    119                 buf=append(buf,zeros((2,windowsize-buf.shape[1])),1)
    120             buf=buf*window
    121     
    122             #get the amplitudes of the frequency components and discard the phases
    123             freqs=abs(fft.rfft(buf))
    124 
    125             #scale down the spectrum to detect onsets
    126             freqs_len=freqs.shape[1]
    127             if num_bins_scaled_freq<freqs_len:
    128                 freqs_len_div=freqs_len//num_bins_scaled_freq
    129                 new_freqs_len=freqs_len_div*num_bins_scaled_freq
    130                 freqs_scaled=mean(mean(freqs,0)[:new_freqs_len].reshape([num_bins_scaled_freq,freqs_len_div]),1)
    131             else:
    132                 freqs_scaled=zeros(num_bins_scaled_freq)
    133 
    134 
    135             #process onsets
    136             m=2.0*mean(freqs_scaled-old_freqs_scaled)/(mean(abs(old_freqs_scaled))+1e-3)
    137             if m<0.0:
    138                 m=0.0
    139             if m>1.0:
    140                 m=1.0
    141             if plot_onsets:
    142                 onsets.append(m)
    143             if m>onset_level:
    144                 displace_tick=1.0
    145                 extra_onset_time_credit+=1.0
    146 
    147         cfreqs=(freqs*displace_tick)+(old_freqs*(1.0-displace_tick))
    148 
    149         #randomize the phases by multiplication with a random complex number with modulus=1
    150         ph=random.uniform(0,2*pi,(nchannels,cfreqs.shape[1]))*1j
    151         cfreqs=cfreqs*exp(ph)
    152 
    153         #do the inverse FFT 
    154         buf=fft.irfft(cfreqs)
    155 
    156         #window again the output buffer
    157         buf*=window
    158 
    159         #overlap-add the output
    160         output=buf[:,0:half_windowsize]+old_windowed_buf[:,half_windowsize:windowsize]
    161         old_windowed_buf=buf
    162 
    163         #remove the resulted amplitude modulation
    164         output*=hinv_buf
    165         
    166         #clamp the values to -1..1 
    167         output[output>1.0]=1.0
    168         output[output<-1.0]=-1.0
    169 
    170         #write the output to wav file
    171         outfile.writeframes(int16(output.ravel(1)*32767.0).tostring())
    172 
    173         if get_next_buf:
    174             start_pos+=displace_pos
    175 
    176         get_next_buf=False
    177 
    178         if start_pos>=nsamples:
    179             print ("100 %")
    180             break
    181         sys.stdout.write ("%d %% \r" % int(100.0*start_pos/nsamples))
    182         sys.stdout.flush()
    183 
    184         
    185         if extra_onset_time_credit<=0.0:
    186             displace_tick+=displace_tick_increase
    187         else:
    188             credit_get=0.5*displace_tick_increase #this must be less than displace_tick_increase
    189             extra_onset_time_credit-=credit_get
    190             if extra_onset_time_credit<0:
    191                 extra_onset_time_credit=0
    192             displace_tick+=displace_tick_increase-credit_get
    193 
    194         if displace_tick>=1.0:
    195             displace_tick=displace_tick % 1.0
    196             get_next_buf=True
    197 
    198     outfile.close()
    199     
    200     if plot_onsets:
    201         plt.plot(onsets)
    202         plt.show()
    203     
    204 
    205 ########################################
    206 print ("Paul's Extreme Sound Stretch (Paulstretch) - Python version 20141220")
    207 print ("new method: using onsets information")
    208 print ("by Nasca Octavian PAUL, Targu Mures, Romania\n")
    209 parser = OptionParser(usage="usage: %prog [options] input_wav output_wav")
    210 parser.add_option("-s", "--stretch", dest="stretch",help="stretch amount (1.0 = no stretch)",type="float",default=8.0)
    211 parser.add_option("-w", "--window_size", dest="window_size",help="window size (seconds)",type="float",default=0.25)
    212 parser.add_option("-t", "--onset", dest="onset",help="onset sensitivity (0.0=max,1.0=min)",type="float",default=10.0)
    213 (options, args) = parser.parse_args()
    214 
    215 
    216 if (len(args)<2) or (options.stretch<=0.0) or (options.window_size<=0.001):
    217     print ("Error in command line parameters. Run this program with --help for help.")
    218     sys.exit(1)
    219 
    220 print ("stretch amount = %g" % options.stretch)
    221 print ("window size = %g seconds" % options.window_size)
    222 print ("onset sensitivity = %g" % options.onset)
    223 (samplerate,smp)=load_wav(args[0])
    224 
    225 paulstretch(samplerate,smp,options.stretch,options.window_size,options.onset,args[1])
    226 
    227 
    228