1import os.path
2import numpy as np
3import itertools
4import Tools
5from scipy import signal
6#from pylab import figure, clf, plot, xlabel, ylabel, xlim, ylim, title, grid, axes, show,semilogx, semilogy
7import math
8
9# Those patterns are used for tests and benchmarks.
10# For tests, there is the need to add tests for saturation
11
12def cartesian(*somelists):
13   r=[]
14   for element in itertools.product(*somelists):
15       r.append(element)
16   return(r)
17
18def writeBenchmarks(config):
19    NBSAMPLES=512 # 512 for stereo
20    NUMSTAGES = 4
21
22    samples=np.random.randn(NBSAMPLES)
23    coefs=np.random.randn(NUMSTAGES*5)
24
25    samples = Tools.normalize(samples)
26    coefs = Tools.normalize(coefs)
27
28
29    # Used for benchmarks
30    config.writeInput(1, samples,"Samples")
31    config.writeInput(1, coefs,"Coefs")
32
33def getCoefs(n,sos,format):
34    if format==15:
35       coefs=np.reshape(np.hstack((np.insert(sos[:,:3],1,0.0,axis=1),-sos[:,4:])),n*6)
36    else:
37       coefs=np.reshape(np.hstack((sos[:,:3],-sos[:,4:])),n*5)
38
39    if format==31:
40        # Postshift must be 2 in the tests
41        coefs = coefs / 4.0
42
43    if format==15:
44        # Postshift must be 2 in the tests
45        coefs = coefs / 4.0
46
47    return(coefs)
48
49def genSos(numTaps):
50    zeros=[]
51    poles=[]
52    for i in range(0,numTaps):
53        phase = np.random.rand()*2.0 * math.pi
54        z = np.exp(1j*phase)
55
56        phase = np.random.rand()*2.0 * math.pi
57        amplitude = np.random.rand()*0.7
58        p = np.exp(1j*phase) * amplitude
59
60        zeros += [z,np.conj(z)]
61        poles += [p,np.conj(p)]
62
63    g = 0.02
64
65    sos = signal.zpk2sos(zeros,poles,g)
66
67    return(sos)
68
69
70def writeTests(config,format):
71    # Write test with fixed and known patterns
72    NB = 100
73    t = np.linspace(0, 1,NB)
74
75    sig = Tools.normalize(np.sin(2*np.pi*5*t)+np.random.randn(len(t)) * 0.2 + 0.4*np.sin(2*np.pi*20*t))
76
77    if format==31:
78       sig = 1.0*sig / (1 << 2)
79
80    #if format==15:
81    #   sig = 1.0*sig / 2.0
82
83    p0 = np.exp(1j*0.05) * 0.98
84    p1 = np.exp(1j*0.25) * 0.9
85    p2 = np.exp(1j*0.45) * 0.97
86
87    z0 = np.exp(1j*0.02)
88    z1 = np.exp(1j*0.65)
89    z2 = np.exp(1j*1.0)
90
91    g = 0.02
92
93    sos = signal.zpk2sos(
94          [z0,np.conj(z0),z1,np.conj(z1),z2,np.conj(z2)]
95         ,[p0, np.conj(p0),p1, np.conj(p1),p2, np.conj(p2)]
96         ,g)
97
98    coefs=getCoefs(3,sos,format)
99
100    res=signal.sosfilt(sos,sig)
101
102    config.writeInput(1, sig,"BiquadInput")
103    config.writeInput(1, res,"BiquadOutput")
104    config.writeInput(1, coefs,"BiquadCoefs")
105
106    #if format==0:
107    #   figure()
108    #   plot(sig)
109    #   figure()
110    #   plot(res)
111    #   show()
112
113    # Now random patterns to test different tail sizes
114    # and number of loops
115
116    numStages = [Tools.loopnb(format,Tools.TAILONLY),
117       Tools.loopnb(format,Tools.BODYONLY),
118       Tools.loopnb(format,Tools.BODYANDTAIL)
119    ]
120
121    blockSize=[Tools.loopnb(format,Tools.TAILONLY),
122       Tools.loopnb(format,Tools.BODYONLY),
123       Tools.loopnb(format,Tools.BODYANDTAIL)
124    ]
125
126    allConfigs = cartesian(numStages, blockSize)
127
128    allconf=[]
129    allcoefs=[]
130    allsamples=[]
131    allStereo=[]
132    alloutputs=[]
133    allStereoOutputs=[]
134
135    for (n,b) in allConfigs:
136        samples=np.random.randn(b)
137        samples = Tools.normalize(samples)
138
139        samplesB=np.random.randn(b)
140        samplesB = Tools.normalize(samplesB)
141
142        stereo = np.empty((samples.size + samplesB.size,), dtype=samples.dtype)
143        stereo[0::2] = samples
144        stereo[1::2] = samplesB
145
146        sos = genSos(n)
147        coefs=getCoefs(n,sos,format)
148
149        output=signal.sosfilt(sos,samples)
150        outputB=signal.sosfilt(sos,samplesB)
151
152        stereoOutput = np.empty((output.size + outputB.size,), dtype=output.dtype)
153        stereoOutput[0::2] = output
154        stereoOutput[1::2] = outputB
155
156        allStereoOutputs += list(stereoOutput)
157        alloutputs += list(output)
158        allconf += [n,b]
159        allcoefs += list(coefs)
160        allsamples += list(samples)
161        allStereo += list(stereo)
162
163
164    config.writeReferenceS16(2,allconf,"AllBiquadConfigs")
165    config.writeInput(2,allsamples,"AllBiquadInputs")
166    config.writeInput(2,allcoefs,"AllBiquadCoefs")
167    config.writeReference(2,alloutputs,"AllBiquadRefs")
168    # Stereo version only for floats
169    if format==0 or format==16:
170        config.writeInput(2,allStereo,"AllBiquadStereoInputs")
171        config.writeReference(2,allStereoOutputs,"AllBiquadStereoRefs")
172
173
174
175def generatePatterns():
176    PATTERNDIR = os.path.join("Patterns","DSP","Filtering","BIQUAD","BIQUAD")
177    PARAMDIR = os.path.join("Parameters","DSP","Filtering","BIQUAD","BIQUAD")
178
179    configf64=Tools.Config(PATTERNDIR,PARAMDIR,"f64")
180    configf32=Tools.Config(PATTERNDIR,PARAMDIR,"f32")
181    configf16=Tools.Config(PATTERNDIR,PARAMDIR,"f16")
182    configq31=Tools.Config(PATTERNDIR,PARAMDIR,"q31")
183    configq15=Tools.Config(PATTERNDIR,PARAMDIR,"q15")
184    #configq7=Tools.Config(PATTERNDIR,PARAMDIR,"q7")
185
186
187
188    writeBenchmarks(configf32)
189    writeBenchmarks(configf16)
190    writeBenchmarks(configq31)
191    writeBenchmarks(configq15)
192    writeBenchmarks(configf64)
193
194    writeTests(configf32,0)
195    writeTests(configf16,16)
196    writeTests(configq31,31)
197    writeTests(configq15,15)
198    writeTests(configf64,64)
199
200    #writeTests(configq7)
201
202if __name__ == '__main__':
203  generatePatterns()
204
205