xref: /btstack/test/sbc/sbc.py (revision f08a674b88cd54d27621208e04242ed1acac4ed1)
1#!/usr/bin/env python
2import numpy as np
3import wave
4import struct
5import sys
6
7# channel mode
8MONO = 0
9DUAL_CHANNEL = 1
10STEREO = 2
11JOINT_STEREO = 3
12channel_modes = ["MONO", "DUAL CHANNEL", "STEREO", "JOINT STEREO"]
13
14# allocation method
15LOUDNESS = 0
16SNR = 1
17allocation_methods = ["LOUDNESS", "SNR"]
18
19sampling_frequencies = [16000, 32000, 44100, 48000]
20nr_blocks = [4, 8, 12, 16]
21nr_subbands = [4, 8]
22
23
24def allocation_method_to_str(allocation_method):
25    global allocation_methods
26    return allocation_methods[allocation_method]
27
28def channel_mode_to_str(channel_mode):
29    global channel_modes
30    return channel_modes[channel_mode]
31
32def sampling_frequency_to_str(sampling_frequency):
33    global sampling_frequencies
34    return sampling_frequencies[sampling_frequency]
35
36def sampling_frequency_index(sampling_frequency):
37    global sampling_frequencies
38    for index, value in enumerate(sampling_frequencies):
39        if value == sampling_frequency:
40            return index
41    return -1
42
43Proto_4_40 = [
44    0.00000000E+00, 5.36548976E-04, 1.49188357E-03, 2.73370904E-03,
45    3.83720193E-03, 3.89205149E-03, 1.86581691E-03, -3.06012286E-03,
46    1.09137620E-02, 2.04385087E-02, 2.88757392E-02, 3.21939290E-02,
47    2.58767811E-02, 6.13245186E-03, -2.88217274E-02, -7.76463494E-02,
48    1.35593274E-01, 1.94987841E-01, 2.46636662E-01, 2.81828203E-01,
49    2.94315332E-01, 2.81828203E-01, 2.46636662E-01, 1.94987841E-01,
50    -1.35593274E-01, -7.76463494E-02, -2.88217274E-02, 6.13245186E-03,
51    2.58767811E-02, 3.21939290E-02, 2.88757392E-02, 2.04385087E-02,
52    -1.09137620E-02, -3.06012286E-03, 1.86581691E-03, 3.89205149E-03,
53    3.83720193E-03, 2.73370904E-03, 1.49188357E-03, 5.36548976E-04
54]
55
56Proto_8_80 = [
57    0.00000000E+00, 1.56575398E-04, 3.43256425E-04, 5.54620202E-04,
58    8.23919506E-04, 1.13992507E-03, 1.47640169E-03, 1.78371725E-03,
59    2.01182542E-03, 2.10371989E-03, 1.99454554E-03, 1.61656283E-03,
60    9.02154502E-04, -1.78805361E-04, -1.64973098E-03, -3.49717454E-03,
61    5.65949473E-03, 8.02941163E-03, 1.04584443E-02, 1.27472335E-02,
62    1.46525263E-02, 1.59045603E-02, 1.62208471E-02, 1.53184106E-02,
63    1.29371806E-02, 8.85757540E-03, 2.92408442E-03, -4.91578024E-03,
64    -1.46404076E-02, -2.61098752E-02, -3.90751381E-02, -5.31873032E-02,
65    6.79989431E-02, 8.29847578E-02, 9.75753918E-02, 1.11196689E-01,
66    1.23264548E-01, 1.33264415E-01, 1.40753505E-01, 1.45389847E-01,
67    1.46955068E-01, 1.45389847E-01, 1.40753505E-01, 1.33264415E-01,
68    1.23264548E-01, 1.11196689E-01, 9.75753918E-02, 8.29847578E-02,
69    -6.79989431E-02, -5.31873032E-02, -3.90751381E-02, -2.61098752E-02,
70    -1.46404076E-02, -4.91578024E-03, 2.92408442E-03, 8.85757540E-03,
71    1.29371806E-02, 1.53184106E-02, 1.62208471E-02, 1.59045603E-02,
72    1.46525263E-02, 1.27472335E-02, 1.04584443E-02, 8.02941163E-03,
73    -5.65949473E-03, -3.49717454E-03, -1.64973098E-03, -1.78805361E-04,
74    9.02154502E-04, 1.61656283E-03, 1.99454554E-03, 2.10371989E-03,
75    2.01182542E-03, 1.78371725E-03, 1.47640169E-03, 1.13992507E-03,
76    8.23919506E-04, 5.54620202E-04, 3.43256425E-04, 1.56575398E-04
77]
78
79crc_table = [
80    0x00, 0x1D, 0x3A, 0x27, 0x74, 0x69, 0x4E, 0x53,
81    0xE8, 0xF5, 0xD2, 0xCF, 0x9C, 0x81, 0xA6, 0xBB,
82    0xCD, 0xD0, 0xF7, 0xEA, 0xB9, 0xA4, 0x83, 0x9E,
83    0x25, 0x38, 0x1F, 0x02, 0x51, 0x4C, 0x6B, 0x76,
84    0x87, 0x9A, 0xBD, 0xA0, 0xF3, 0xEE, 0xC9, 0xD4,
85    0x6F, 0x72, 0x55, 0x48, 0x1B, 0x06, 0x21, 0x3C,
86    0x4A, 0x57, 0x70, 0x6D, 0x3E, 0x23, 0x04, 0x19,
87    0xA2, 0xBF, 0x98, 0x85, 0xD6, 0xCB, 0xEC, 0xF1,
88    0x13, 0x0E, 0x29, 0x34, 0x67, 0x7A, 0x5D, 0x40,
89    0xFB, 0xE6, 0xC1, 0xDC, 0x8F, 0x92, 0xB5, 0xA8,
90    0xDE, 0xC3, 0xE4, 0xF9, 0xAA, 0xB7, 0x90, 0x8D,
91    0x36, 0x2B, 0x0C, 0x11, 0x42, 0x5F, 0x78, 0x65,
92    0x94, 0x89, 0xAE, 0xB3, 0xE0, 0xFD, 0xDA, 0xC7,
93    0x7C, 0x61, 0x46, 0x5B, 0x08, 0x15, 0x32, 0x2F,
94    0x59, 0x44, 0x63, 0x7E, 0x2D, 0x30, 0x17, 0x0A,
95    0xB1, 0xAC, 0x8B, 0x96, 0xC5, 0xD8, 0xFF, 0xE2,
96    0x26, 0x3B, 0x1C, 0x01, 0x52, 0x4F, 0x68, 0x75,
97    0xCE, 0xD3, 0xF4, 0xE9, 0xBA, 0xA7, 0x80, 0x9D,
98    0xEB, 0xF6, 0xD1, 0xCC, 0x9F, 0x82, 0xA5, 0xB8,
99    0x03, 0x1E, 0x39, 0x24, 0x77, 0x6A, 0x4D, 0x50,
100    0xA1, 0xBC, 0x9B, 0x86, 0xD5, 0xC8, 0xEF, 0xF2,
101    0x49, 0x54, 0x73, 0x6E, 0x3D, 0x20, 0x07, 0x1A,
102    0x6C, 0x71, 0x56, 0x4B, 0x18, 0x05, 0x22, 0x3F,
103    0x84, 0x99, 0xBE, 0xA3, 0xF0, 0xED, 0xCA, 0xD7,
104    0x35, 0x28, 0x0F, 0x12, 0x41, 0x5C, 0x7B, 0x66,
105    0xDD, 0xC0, 0xE7, 0xFA, 0xA9, 0xB4, 0x93, 0x8E,
106    0xF8, 0xE5, 0xC2, 0xDF, 0x8C, 0x91, 0xB6, 0xAB,
107    0x10, 0x0D, 0x2A, 0x37, 0x64, 0x79, 0x5E, 0x43,
108    0xB2, 0xAF, 0x88, 0x95, 0xC6, 0xDB, 0xFC, 0xE1,
109    0x5A, 0x47, 0x60, 0x7D, 0x2E, 0x33, 0x14, 0x09,
110    0x7F, 0x62, 0x45, 0x58, 0x0B, 0x16, 0x31, 0x2C,
111    0x97, 0x8A, 0xAD, 0xB0, 0xE3, 0xFE, 0xD9, 0xC4
112]
113# offset table for 4-subbands
114offset4 = np.array([[ -1, 0, 0, 0 ],
115                    [ -2, 0, 0, 1 ],
116                    [ -2, 0, 0, 1 ],
117                    [ -2, 0, 0, 1 ]
118                    ])
119
120# offset tables for 8-subbands
121offset8 = np.array([[ -2, 0, 0, 0, 0, 0, 0, 1 ],
122                    [ -3, 0, 0, 0, 0, 0, 1, 2 ],
123                    [ -4, 0, 0, 0, 0, 0, 1, 2 ],
124                    [ -4, 0, 0, 0, 0, 0, 1, 2 ]
125                    ])
126
127def calculate_scalefactor(max_subbandsample):
128    x = 0
129    while True:
130        y = 1 << x + 1
131        if y > max_subbandsample:
132            break
133        x += 1
134    return (x,y)
135
136
137def calculate_max_subbandsample(nr_blocks, nr_channels, nr_subbands, sb_sample):
138    max_subbandsample = np.zeros(shape = (nr_channels, nr_subbands))
139
140    for blk in range(nr_blocks):
141        for ch in range(nr_channels):
142            for sb in range(nr_subbands):
143                m = abs(sb_sample[blk][ch][sb])
144                if max_subbandsample[ch][sb] < m:
145                    max_subbandsample[ch][sb] = m
146    return max_subbandsample
147
148def calculate_scalefactors(nr_blocks, nr_channels, nr_subbands, sb_sample):
149    scale_factor =  np.zeros(shape=(nr_channels, nr_subbands), dtype = np.int32)
150    scalefactor =  np.zeros(shape=(nr_channels, nr_subbands), dtype = np.int32)
151
152    # max_subbandsample = calculate_max_subbandsample(nr_blocks, nr_channels, nr_subbands, sb_sample)
153    # for ch in range(nr_channels):
154    #     for sb in range(nr_subbands):
155    #         (scale_factor[ch][sb], scalefactor[ch][sb]) = calculate_scalefactor(max_subbandsample[ch][sb])
156
157    for ch in range(nr_channels):
158        for sb in range(nr_subbands):
159            scale_factor[ch][sb] = 0
160            scalefactor[ch][sb] = 2
161            for blk in range(nr_blocks):
162                while scalefactor[ch][sb] < abs(sb_sample[blk][ch][sb]):
163                    scale_factor[ch][sb]+=1
164                    scalefactor[ch][sb] *= 2
165
166    return scale_factor, scalefactor
167
168def calculate_scalefactors_and_channel_mode(frame):
169    frame.scale_factor, frame.scalefactor = calculate_scalefactors(frame.nr_blocks, frame.nr_channels, frame.nr_subbands, frame.sb_sample)
170    #print "calculate_scalefactors_and_channel_mode1 ", frame.scale_factor
171
172    if frame.nr_channels == 1:
173        frame.channel_mode = MONO
174    else:
175        sb_sample1 = np.zeros(shape = (frame.nr_blocks,2,frame.nr_subbands), dtype = np.uint16)
176
177        for blk in range(frame.nr_blocks):
178            for sb in range(frame.nr_subbands):
179                sb_sample1[blk][0][sb] = frame.sb_sample[blk][0][sb] + frame.sb_sample[blk][1][sb]
180                sb_sample1[blk][1][sb] = frame.sb_sample[blk][0][sb] - frame.sb_sample[blk][1][sb]
181
182        scale_factor, scalefactor = calculate_scalefactors(frame.nr_blocks, frame.nr_channels, frame.nr_subbands, sb_sample1)
183        #print "calculate_scalefactors_and_channel_mode 2", scale_factor
184        sumb = 0
185        suma = 0
186        for sb in range(frame.nr_subbands):
187            suma += frame.scale_factor[0][sb] + frame.scale_factor[1][sb]
188            sumb += scale_factor[0][sb] + scale_factor[1][sb]
189
190        #print "calculate_scalefactors_and_channel_mode 3", suma, sumb
191        if suma > sumb:
192            frame.channel_mode = JOINT_STEREO
193        else:
194            frame.channel_mode = STEREO
195
196
197class SBCFrame:
198    syncword = 0
199    sampling_frequency = 0
200    nr_blocks = 0
201    channel_mode = 0
202    nr_channels = 0
203    allocation_method = 0
204    nr_subbands = 0
205    bitpool = 0
206    crc_check = 0
207    # pro subband - 1
208    join = np.zeros(8, dtype = np.uint8)
209    scale_factor =  np.zeros(shape=(2, 8), dtype = np.int32)
210    scalefactor =  np.zeros(shape=(2, 8), dtype = np.int32)
211    audio_sample = np.zeros(shape = (16,2,8), dtype = np.uint16)
212    sb_sample = np.zeros(shape = (16,2,8), dtype = np.uint16)
213    X = np.zeros(8, dtype = np.int16)
214    EX = np.zeros(8)
215    pcm = np.zeros(shape=(2, 8*16), dtype = np.int16)
216    bits    = np.zeros(shape=(2, 8))
217    levels = np.zeros(shape=(2, 8), dtype = np.int32)
218
219
220    def __init__(self, nr_blocks=16, nr_subbands=4, nr_channels=1, bitpool=31, sampling_frequency=44100, allocation_method = 0):
221        self.nr_blocks = nr_blocks
222        self.nr_subbands = nr_subbands
223        self.nr_channels = nr_channels
224        self.sampling_frequency = sampling_frequency_index(sampling_frequency)
225        self.bitpool = bitpool
226        self.allocation_method = allocation_method
227        self.scale_factor = np.zeros(shape=(nr_channels, nr_subbands), dtype = np.int32)
228        self.scalefactor = np.zeros(shape=(nr_channels, nr_subbands), dtype = np.int32)
229        self.audio_sample = np.zeros(shape=(nr_blocks, nr_channels, nr_subbands), dtype = np.uint16)
230        self.sb_sample = np.zeros(shape=(nr_blocks, nr_channels, nr_subbands), dtype = np.uint16)
231        self.levels = np.zeros(shape=(nr_channels, nr_subbands), dtype = np.int32)
232        self.EX = np.zeros(nr_subbands)
233        return
234
235    def dump_audio_samples(self, blk, ch):
236        print self.audio_sample[blk][ch]
237
238    def dump_subband_samples(self, blk, ch):
239        print self.sb_sample[blk][ch]
240
241    def dump_state(self):
242        res =  "SBCFrameHeader state:"
243        res += "\n - nr channels %d" % self.nr_channels
244        res += "\n - nr blocks %d" % self.nr_blocks
245        res += "\n - nr subbands %d" % self.nr_subbands
246        res += "\n - scale factors: %s" % self.scale_factor
247        res += "\n - levels: %s" % self.levels
248        res += "\n - join: %s" % self.join
249        res += "\n - bits: %s" % self.bits
250        print res
251
252    def __str__(self):
253        res =  "SBCFrameHeader:"
254        res += "\n - syncword %d" % self.syncword
255        res += "\n - sampling frequency %d Hz" % sampling_frequency_to_str(self.sampling_frequency)
256
257        res += "\n - nr channels %d" % self.nr_channels
258        res += "\n - nr blocks %d" % self.nr_blocks
259        res += "\n - nr subbands %d" % self.nr_subbands
260
261        res += "\n - channel mode %s" % channel_mode_to_str(self.channel_mode)
262        res += "\n - allocation method %s" % allocation_method_to_str(self.allocation_method)
263
264        res += "\n - bitpool %d" % self.bitpool
265        res += "\n - crc check %d" % self.crc_check
266        return res
267
268
269def sbc_bit_allocation_stereo_joint(frame):
270    bitneed = np.zeros(shape=(frame.nr_channels, frame.nr_subbands), dtype = np.int32)
271    bits    = np.zeros(shape=(frame.nr_channels, frame.nr_subbands), dtype = np.int32)
272
273    loudness = 0
274
275    if frame.allocation_method == SNR:
276        for ch in range(frame.nr_channels):
277            for sb in range(frame.nr_subbands):
278                bitneed[ch][sb] = frame.scale_factor[ch][sb]
279    else:
280        for ch in range(frame.nr_channels):
281            for sb in range(frame.nr_subbands):
282                if frame.scale_factor[ch][sb] == 0:
283                    bitneed[ch][sb] = -5
284                else:
285                    if frame.nr_subbands == 4:
286                        loudness = scale_factor[ch][sb] - offset4[frame.sampling_frequency][sb]
287                    else:
288                        loudness = frame.scale_factor[ch][sb] - offset8[frame.sampling_frequency][sb]
289
290                    if loudness > 0:
291                        bitneed[ch][sb] = loudness/2
292                    else:
293                        bitneed[ch][sb] = loudness
294
295    # search the maximum bitneed index
296    max_bitneed = 0
297    for ch in range(frame.nr_channels):
298        for sb in range(frame.nr_subbands):
299            if bitneed[ch][sb] > max_bitneed:
300                max_bitneed = bitneed[ch][sb]
301
302    # calculate how many bitslices fit into the bitpool
303    bitcount = 0
304    slicecount = 0
305    bitslice = max_bitneed + 1 #/* init just above the largest sf */
306
307    while True:
308        bitslice -= 1
309        bitcount += slicecount
310        slicecount = 0
311        for ch in range(frame.nr_channels):
312            for sb in range(frame.nr_subbands):
313                if (bitneed[ch][sb] > bitslice+1) and (bitneed[ch][sb] < bitslice+16):
314                    slicecount += 1
315                elif bitneed[ch][sb] == bitslice + 1:
316                    slicecount += 2
317        if bitcount + slicecount >= frame.bitpool:
318            break
319
320    if bitcount + slicecount == frame.bitpool:
321        bitcount += slicecount
322        bitslice -= 1
323
324    # bits are distributed until the last bitslice is reached
325    for ch in range(frame.nr_channels):
326        for sb in range(frame.nr_subbands):
327            if bitneed[ch][sb] < bitslice+2 :
328               bits[ch][sb]=0;
329            else:
330                bits[ch][sb] = min(bitneed[ch][sb]-bitslice,16)
331
332    ch = 0
333    sb = 0
334    while bitcount < frame.bitpool and sb < frame.nr_subbands:
335        if bits[ch][sb] >= 2 and bits[ch][sb] < 16:
336            bits[ch][sb] += 1
337            bitcount += 1
338        elif (bitneed[ch][sb] == bitslice+1) and (frame.bitpool > bitcount+1):
339            bits[ch][sb] = 2
340            bitcount += 2
341        if ch == 1:
342            ch = 0
343            sb += 1
344        else:
345            ch = 1
346
347
348    ch = 0
349    sb = 0
350    while bitcount < frame.bitpool and sb < frame.nr_subbands:
351        if bits[ch][sb] < 16:
352            bits[ch][sb]+=1
353            bitcount+=1
354        if ch == 1:
355            ch = 0
356            sb += 1
357        else:
358            ch = 1
359
360    if bits.sum() != frame.bitpool:
361        print "bit allocation failed, bitpool %d, allocated %d" % (bits.sum() , frame.bitpool)
362        exit(1)
363    return bits
364
365
366def sbc_bit_allocation_mono_dual(frame):
367    #print "Bit allocation for mono/dual channel"
368    bitneed = np.zeros(shape=(frame.nr_channels, frame.nr_subbands), dtype = np.int32)
369    bits    = np.zeros(shape=(frame.nr_channels, frame.nr_subbands), dtype = np.int32)
370    loudness = 0
371
372    for ch in range(frame.nr_channels):
373        # bitneed values are derived from the scale factors
374        if frame.allocation_method == SNR:
375            for sb in range(frame.nr_subbands):
376                bitneed[ch][sb] = frame.scale_factor[ch][sb]
377        else:
378            for sb in range(frame.nr_subbands):
379                if frame.scale_factor[ch][sb] == 0:
380                    bitneed[ch][sb] = -5
381                else:
382                    if frame.nr_subbands == 4:
383                        loudness = frame.scale_factor[ch][sb] - offset4[frame.sampling_frequency][sb]
384                    else:
385                        loudness = frame.scale_factor[ch][sb] - offset8[frame.sampling_frequency][sb]
386                    if loudness > 0:
387                        bitneed[ch][sb] = loudness/2
388                    else:
389                        bitneed[ch][sb] = loudness
390
391        # search the maximum bitneed index
392        max_bitneed = 0
393        for sb in range(frame.nr_subbands):
394            if bitneed[ch][sb] > max_bitneed:
395                max_bitneed = bitneed[ch][sb]
396
397        # calculate how many bitslices fit into the bitpool
398        bitcount = 0
399        slicecount = 0
400        bitslice = max_bitneed + 1
401
402        while True:
403            bitslice = bitslice - 1
404            bitcount = bitcount + slicecount
405            slicecount = 0
406            for sb in range(frame.nr_subbands):
407                if (bitneed[ch][sb] > bitslice+1) and (bitneed[ch][sb] < bitslice+16):
408                    slicecount = slicecount + 1
409                elif bitneed[ch][sb] == bitslice + 1:
410                    slicecount = slicecount + 2
411            if bitcount + slicecount >= frame.bitpool:
412                break
413
414        if bitcount + slicecount == frame.bitpool:
415            bitcount = bitcount + slicecount
416            bitslice = bitslice - 1
417
418        for sb in range(frame.nr_subbands):
419            if bitneed[ch][sb] < bitslice+2 :
420               bits[ch][sb]=0;
421            else:
422                bits[ch][sb] = min(bitneed[ch][sb]-bitslice,16)
423
424        sb = 0
425        while bitcount < frame.bitpool and sb < frame.nr_subbands:
426            if bits[ch][sb] >= 2 and bits[ch][sb] < 16:
427                   bits[ch][sb] = bits[ch][sb] + 1
428                   bitcount = bitcount + 1
429
430            elif (bitneed[ch][sb] == bitslice+1) and (frame.bitpool > bitcount+1):
431                bits[ch][sb] = 2
432                bitcount += 2
433
434            sb = sb + 1
435
436
437        sb = 0
438        while bitcount < frame.bitpool and sb < frame.nr_subbands:
439            if bits[ch][sb] < 16:
440                bits[ch][sb] = bits[ch][sb] + 1
441                bitcount = bitcount + 1
442            sb = sb + 1
443
444    return bits
445
446def sbc_bit_allocation(frame):
447    if frame.channel_mode == MONO or frame.channel_mode == DUAL_CHANNEL:
448        return sbc_bit_allocation_mono_dual(frame)
449    elif frame.channel_mode == STEREO or frame.channel_mode == JOINT_STEREO:
450        return sbc_bit_allocation_stereo_joint(frame)
451    else:
452        print "Wrong channel mode ", frame.channel_mode
453        return -1
454
455def sbc_sampling_frequency_index(sample_rate):
456    sbc_sampling_frequency_index = 0
457    for i in range(len(sampling_frequency)):
458        if sample_rate == sampling_frequency[i]:
459            sbc_sampling_frequency_index = i
460            break
461    return sbc_sampling_frequency_index
462
463
464def sbc_crc8(data, data_len):
465    crc = 0x0f
466    j = 0
467    for i in range(data_len / 8):
468        crc = crc_table[crc ^ data[i]]
469        j = i + 1
470
471    bits_left = data_len%8
472    if bits_left:
473        octet = data[j]
474        for i in range(data_len%8):
475            bit = ((octet ^ crc) & 0x80) >> 7
476            if bit:
477                bit = 0x1d
478            crc = ((crc & 0x7f) << 1) ^ bit
479            octet = octet << 1
480    return crc
481
482
483bitstream = None
484bitstream_index = -1
485bitstream_bits_available = 0
486
487def init_bitstream():
488    global bitstream, bitstream_bits_available, bitstream_index
489    bitstream = []
490    bitstream_index = -1
491    bitstream_bits_available = 0
492
493def add_bit(bit):
494    global bitstream, bitstream_bits_available, bitstream_index
495    if bitstream_bits_available == 0:
496        bitstream.append(0)
497        bitstream_bits_available = 8
498        bitstream_index += 1
499
500    bitstream[bitstream_index] |= bit << (bitstream_bits_available - 1)
501    bitstream_bits_available -= 1
502
503
504def add_bits(bits, len):
505    global bitstream, bitstream_bits_available
506    for i in range(len):
507        add_bit((bits >> (len-1-i)) & 1)
508
509ibuffer = None
510ibuffer_count = 0
511
512def get_bit(fin):
513    global ibuffer, ibuffer_count
514    if ibuffer_count == 0:
515        ibuffer = ord(fin.read(1))
516        ibuffer_count = 8
517
518    bit = (ibuffer >> 7) & 1
519    ibuffer = ibuffer << 1
520    ibuffer_count = ibuffer_count - 1
521    return bit
522
523def drop_remaining_bits():
524    global ibuffer_count
525    #print "dropping %d bits" % ibuffer_count
526    ibuffer_count = 0
527
528def get_bits(fin, bit_count):
529    bits = 0
530    for i in range(bit_count):
531        bits = (bits << 1) | get_bit(fin)
532    # print "get bits: %d -> %02x" %(bit_count, bits)
533    return bits
534
535
536def calculate_crc(frame):
537    global bitstream, bitstream_bits_available, bitstream_index
538    init_bitstream()
539
540    add_bits(frame.sampling_frequency, 2)
541    add_bits(frame.nr_blocks/4-1, 2)
542    add_bits(frame.channel_mode, 2)
543    add_bits(frame.allocation_method, 1)
544    add_bits(frame.nr_subbands/4-1, 1)
545    add_bits(frame.bitpool, 8)
546
547    if frame.channel_mode == JOINT_STEREO:
548        for sb in range(frame.nr_subbands):
549            add_bits(frame.join[sb],1)
550
551    for ch in range(frame.nr_channels):
552        for sb in range(frame.nr_subbands):
553            add_bits(frame.scale_factor[ch][sb], 4)
554
555    bitstream_len = (bitstream_index + 1) * 8
556    if bitstream_bits_available:
557        bitstream_len += (8-bitstream_bits_available)
558    return sbc_crc8(bitstream, bitstream_len)
559
560
561
562def frame_to_bitstream(frame):
563    global bitstream, bitstream_bits_available, bitstream_index
564    init_bitstream()
565
566    add_bits(frame.syncword, 8)
567    add_bits(frame.sampling_frequency, 2)
568    add_bits(frame.nr_blocks/4-1, 2)
569    add_bits(frame.channel_mode, 2)
570    add_bits(frame.allocation_method, 1)
571    add_bits(frame.nr_subbands/4-1, 1)
572    add_bits(frame.bitpool, 8)
573    add_bits(frame.crc_check, 8)
574
575    if frame.channel_mode == JOINT_STEREO:
576        for sb in range(frame.nr_subbands-1):
577            add_bits(frame.join[sb],1)
578        add_bits(0,1)
579
580    for ch in range(frame.nr_channels):
581        for sb in range(frame.nr_subbands):
582            add_bits(frame.scale_factor[ch][sb], 4)
583
584    for blk in range(frame.nr_blocks):
585        for ch in range(frame.nr_channels):
586            for sb in range(frame.nr_subbands):
587                add_bits(frame.audio_sample[blk][ch][sb], frame.bits[ch][sb])
588
589    bitstream_bits_available = 0
590    return bitstream
591
592def mse(a,b):
593    count = 1
594    for i in a.shape:
595        count *= i
596    delta = a - b
597    sqr = delta ** 2
598    res = sqr.sum()*1.0/count
599    # res = ((a - b) ** 2).mean()
600    return res
601