File: examples-hears_time_varying_filter1.txt

package info (click to toggle)
brian 1.4.3-1
  • links: PTS, VCS
  • area: main
  • in suites: sid, stretch
  • size: 23,436 kB
  • sloc: python: 68,707; cpp: 29,040; ansic: 5,182; sh: 111; makefile: 61
file content (122 lines) | stat: -rw-r--r-- 4,497 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
.. currentmodule:: brian

.. index::
   pair: example usage; cos
   pair: example usage; LinearFilterbank
   pair: example usage; log
   pair: example usage; FunctionFilterbank
   pair: example usage; show
   pair: example usage; flipud
   pair: example usage; imshow
   pair: example usage; specgram
   pair: example usage; sqrt
   pair: example usage; sinh
   pair: example usage; whitenoise
   pair: example usage; ControlFilterbank
   pair: example usage; figure
   pair: example usage; squeeze
   pair: example usage; sin
   pair: example usage; arcsinh

.. _example-hears_time_varying_filter1:

Example: time_varying_filter1 (hears)
=====================================

This example implements a band pass filter whose center frequency is modulated
by an Ornstein-Uhlenbeck. The white noise term used for this process is output
by a FunctionFilterbank. The bandpass filter coefficients update is an example
of how to use a :class:`~brian.hears.ControlFilterbank`. The bandpass filter is
a basic biquadratic filter for which the Q factor and the center frequency must
be given. The input is a white noise.

::

    
    from brian import *
    from brian.hears import *
    
    samplerate = 20*kHz
    SoundDuration = 300*ms
    sound = whitenoise(SoundDuration, samplerate).ramp() 
    
    #number of frequency channel (here it must be one as a spectrogram of the
    #output is plotted)
    nchannels = 1   
    
    fc_init = 5000*Hz   #initial center frequency of the band pass filter
    Q = 5               #quality factor of the band pass filter
    update_interval = 4 # the filter coefficients are updated every 4 samples
    
    #parameters of the Ornstein-Uhlenbeck process
    s_i = 1200*Hz
    tau_i = 100*ms      
    mu_i = fc_init/tau_i
    sigma_i = sqrt(2)*s_i/sqrt(tau_i)
    deltaT = defaultclock.dt
    
    #this function  is used in a FunctionFilterbank. It outputs a noise term that
    #will be later used by the controler to update the center frequency
    noise = lambda x: mu_i*deltaT+sigma_i*randn(1)*sqrt(deltaT)
    noise_generator = FunctionFilterbank(sound, noise)
    
    #this class will take as input the output of the noise generator and as target
    #the bandpass filter center frequency
    class CoeffController(object):
        def __init__(self, target):
            self.target = target
            self.deltaT = 1./samplerate
            self.BW = 2*arcsinh(1./2/Q)*1.44269
            self.fc = fc_init
            
        def __call__(self, input):
            #the control variables are taken as the last of the buffer
            noise_term = input[-1,:]
            #update the center frequency by updateing the OU process
            self.fc = self.fc-self.fc/tau_i*self.deltaT+noise_term
    
            w0 = 2*pi*self.fc/samplerate
            #update the coefficient of the biquadratic filterbank
            alpha = sin(w0)*sinh(log(2)/2*self.BW*w0/sin(w0))
            self.target.filt_b[:, 0, 0] = sin(w0)/2
            self.target.filt_b[:, 1, 0] = 0
            self.target.filt_b[:, 2, 0] = -sin(w0)/2
        
            self.target.filt_a[:, 0, 0] = 1+alpha
            self.target.filt_a[:, 1, 0] = -2*cos(w0)
            self.target.filt_a[:, 2, 0] = 1-alpha
            
    # In the present example the time varying filter is a LinearFilterbank therefore
    #we must initialise the filter coefficients; the one used for the first buffer computation
    w0 = 2*pi*fc_init/samplerate
    BW = 2*arcsinh(1./2/Q)*1.44269
    alpha = sin(w0)*sinh(log(2)/2*BW*w0/sin(w0))
    
    filt_b = zeros((nchannels, 3, 1))
    filt_a = zeros((nchannels, 3, 1))
    filt_b[:, 0, 0] = sin(w0)/2
    filt_b[:, 1, 0] = 0
    filt_b[:, 2, 0] = -sin(w0)/2
    filt_a[:, 0, 0] = 1+alpha
    filt_a[:, 1, 0] = -2*cos(w0)
    filt_a[:, 2, 0] = 1-alpha
    
    #the filter which will have time varying coefficients
    bandpass_filter = LinearFilterbank(sound, filt_b, filt_a)
    #the updater
    updater = CoeffController(bandpass_filter)
    
    #the controller. Remember it must be the last of the chain
    control = ControlFilterbank(bandpass_filter, noise_generator, bandpass_filter,
                                updater, update_interval)          
    
    time_varying_filter_mon = control.process()
    
    figure(1)
    pxx, freqs, bins, im = specgram(squeeze(time_varying_filter_mon),
                                    NFFT=256, Fs=samplerate, noverlap=240) 
    imshow(flipud(pxx), aspect='auto')
    
    show()