1# Lint as: python2, python3 2# Copyright 2016 The Chromium OS Authors. All rights reserved. 3# Use of this source code is governed by a BSD-style license that can be 4# found in the LICENSE file. 5 6"""This module provides utilities to detect some artifacts and measure the 7 quality of audio.""" 8 9from __future__ import absolute_import 10from __future__ import division 11from __future__ import print_function 12import logging 13import math 14import numpy 15from six.moves import range 16 17# Normal autotest environment. 18try: 19 import common 20 from autotest_lib.client.cros.audio import audio_analysis 21# Standalone execution without autotest environment. 22except ImportError: 23 import audio_analysis 24 25 26# The input signal should be one sine wave with fixed frequency which 27# can have silence before and/or after sine wave. 28# For example: 29# silence sine wave silence 30# -----------|VVVVVVVVVVVVV|----------- 31# (a) (b) (c) 32# This module detects these artifacts: 33# 1. Detect noise in (a) and (c). 34# 2. Detect delay in (b). 35# 3. Detect burst in (b). 36# Assume the transitions between (a)(b) and (b)(c) are smooth and 37# amplitude increases/decreases linearly. 38# This module will detect artifacts in the sine wave. 39# This module also estimates the equivalent noise level by teager operator. 40# This module also detects volume changes in the sine wave. However, volume 41# changes may be affected by delay or burst. 42# Some artifacts may cause each other. 43 44# In this module, amplitude and frequency are derived from Hilbert transform. 45# Both amplitude and frequency are a function of time. 46 47# To detect each artifact, each point will be compared with 48# average amplitude of its block. The block size will be 1.5 ms. 49# Using average amplitude can mitigate the error caused by 50# Hilbert transform and noise. 51# In some case, for more accuracy, the block size may be modified 52# to other values. 53DEFAULT_BLOCK_SIZE_SECS = 0.0015 54 55# If the difference between average frequency of this block and 56# dominant frequency of full signal is less than 0.5 times of 57# dominant frequency, this block is considered to be within the 58# sine wave. In most cases, if there is no sine wave(only noise), 59# average frequency will be much greater than 5 times of 60# dominant frequency. 61# Also, for delay during playback, the frequency will be about 0 62# in perfect situation or much greater than 5 times of dominant 63# frequency if it's noised. 64DEFAULT_FREQUENCY_ERROR = 0.5 65 66# If the amplitude of some sample is less than 0.6 times of the 67# average amplitude of its left/right block, it will be considered 68# as a delay during playing. 69DEFAULT_DELAY_AMPLITUDE_THRESHOLD = 0.6 70 71# If the average amplitude of the block before or after playing 72# is more than 0.5 times to the average amplitude of the wave, 73# it will be considered as a noise artifact. 74DEFAULT_NOISE_AMPLITUDE_THRESHOLD = 0.5 75 76# In the sine wave, if the amplitude is more than 1.4 times of 77# its left side and its right side, it will be considered as 78# a burst. 79DEFAULT_BURST_AMPLITUDE_THRESHOLD = 1.4 80 81# When detecting burst, if the amplitude is lower than 0.5 times 82# average amplitude, we ignore it. 83DEFAULT_BURST_TOO_SMALL = 0.5 84 85# For a signal which is the combination of sine wave with fixed frequency f and 86# amplitude 1 and standard noise with amplitude k, the average teager value is 87# nearly linear to the noise level k. 88# Given frequency f, we simulate a sine wave with default noise level and 89# calculate its average teager value. Then, we can estimate the equivalent 90# noise level of input signal by the average teager value of input signal. 91DEFAULT_STANDARD_NOISE = 0.005 92 93# For delay, burst, volume increasing/decreasing, if two delay( 94# burst, volume increasing/decreasing) happen within 95# DEFAULT_SAME_EVENT_SECS seconds, we consider they are the 96# same event. 97DEFAULT_SAME_EVENT_SECS = 0.001 98 99# When detecting increasing/decreasing volume of signal, if the amplitude 100# is lower than 0.1 times average amplitude, we ignore it. 101DEFAULT_VOLUME_CHANGE_TOO_SMALL = 0.1 102 103# If average amplitude of right block is less/more than average 104# amplitude of left block times DEFAULT_VOLUME_CHANGE_AMPLITUDE, it will be 105# considered as decreasing/increasing on volume. 106DEFAULT_VOLUME_CHANGE_AMPLITUDE = 0.1 107 108# If the increasing/decreasing volume event is too close to the start or the end 109# of sine wave, we consider its volume change as part of rising/falling phase in 110# the start/end. 111NEAR_START_OR_END_SECS = 0.01 112 113# After applying Hilbert transform, the resulting amplitude and frequency may be 114# extremely large in the start and/or the end part. Thus, we will append zeros 115# before and after the whole wave for 0.1 secs. 116APPEND_ZEROS_SECS = 0.1 117 118# If the noise event is too close to the start or the end of the data, we 119# consider its noise as part of artifacts caused by edge effect of Hilbert 120# transform. 121# For example, originally, the data duration is 10 seconds. 122# We append 0.1 seconds of zeros in the beginning and the end of the data, so 123# the data becomes 10.2 seocnds long. 124# Then, we apply Hilbert transform to 10.2 seconds of data. 125# Near 0.1 seconds and 10.1 seconds, there will be edge effect of Hilbert 126# transform. We do not want these be treated as noise. 127# If NEAR_DATA_START_OR_END_SECS is set to 0.01, then the noise happened 128# at [0, 0.11] and [10.09, 10.1] will be ignored. 129NEAR_DATA_START_OR_END_SECS = 0.01 130 131# If the noise event is too close to the start or the end of the sine wave in 132# the data, we consider its noise as part of artifacts caused by edge effect of 133# Hilbert transform. 134# A |-------------|vvvvvvvvvvvvvvvvvvvvvvv|-------------| 135# B |ooooooooo| d | | d |ooooooooo| 136# 137# A is full signal. It contains a sine wave and silence before and after sine 138# wave. 139# In B, |oooo| shows the parts that we are going to check for noise before/after 140# sine wave. | d | is determined by NEAR_SINE_START_OR_END_SECS. 141NEAR_SINE_START_OR_END_SECS = 0.01 142 143 144class SineWaveNotFound(Exception): 145 """Error when there's no sine wave found in the signal""" 146 pass 147 148 149def hilbert(x): 150 """Hilbert transform copied from scipy. 151 152 More information can be found here: 153 http://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.hilbert.html 154 155 @param x: Real signal data to transform. 156 157 @returns: Analytic signal of x, we can further extract amplitude and 158 frequency from it. 159 160 """ 161 x = numpy.asarray(x) 162 if numpy.iscomplexobj(x): 163 raise ValueError("x must be real.") 164 axis = -1 165 N = x.shape[axis] 166 if N <= 0: 167 raise ValueError("N must be positive.") 168 169 Xf = numpy.fft.fft(x, N, axis=axis) 170 h = numpy.zeros(N) 171 if N % 2 == 0: 172 h[0] = h[N // 2] = 1 173 h[1:N // 2] = 2 174 else: 175 h[0] = 1 176 h[1:(N + 1) // 2] = 2 177 178 if len(x.shape) > 1: 179 ind = [newaxis] * x.ndim 180 ind[axis] = slice(None) 181 h = h[ind] 182 x = numpy.fft.ifft(Xf * h, axis=axis) 183 return x 184 185 186def noised_sine_wave(frequency, rate, noise_level): 187 """Generates a sine wave of 2 second with specified noise level. 188 189 @param frequency: Frequency of sine wave. 190 @param rate: Sampling rate. 191 @param noise_level: Required noise level. 192 193 @returns: A sine wave with specified noise level. 194 195 """ 196 wave = [] 197 for index in range(0, rate * 2): 198 sample = 2.0 * math.pi * frequency * float(index) / float(rate) 199 sine_wave = math.sin(sample) 200 noise = noise_level * numpy.random.standard_normal() 201 wave.append(sine_wave + noise) 202 return wave 203 204 205def average_teager_value(wave, amplitude): 206 """Computes the normalized average teager value. 207 208 After averaging the teager value, we will normalize the value by 209 dividing square of amplitude. 210 211 @param wave: Wave to apply teager operator. 212 @param amplitude: Average amplitude of given wave. 213 214 @returns: Average teager value. 215 216 """ 217 teager_value, length = 0, len(wave) 218 for i in range(1, length - 1): 219 ith_teager_value = abs(wave[i] * wave[i] - wave[i - 1] * wave[i + 1]) 220 ith_teager_value *= max(1, abs(wave[i])) 221 teager_value += ith_teager_value 222 teager_value = (float(teager_value) / length) / (amplitude ** 2) 223 return teager_value 224 225 226def noise_level(amplitude, frequency, rate, teager_value_of_input): 227 """Computes the noise level compared with standard_noise. 228 229 For a signal which is the combination of sine wave with fixed frequency f 230 and amplitude 1 and standard noise with amplitude k, the average teager 231 value is nearly linear to the noise level k. 232 Thus, we can compute the average teager value of a sine wave with 233 standard_noise. Then, we can estimate the noise level of given input. 234 235 @param amplitude: Amplitude of input audio. 236 @param frequency: Dominant frequency of input audio. 237 @param rate: Sampling rate. 238 @param teager_value_of_input: Average teager value of input audio. 239 240 @returns: A float value denotes the audio is equivalent to have 241 how many times of noise compared with its amplitude. 242 For example, 0.02 denotes that the wave has a noise which 243 has standard distribution with standard deviation being 244 0.02 times the amplitude of the wave. 245 246 """ 247 standard_noise = DEFAULT_STANDARD_NOISE 248 249 # Generates the standard sine wave with stdandard_noise level of noise. 250 standard_wave = noised_sine_wave(frequency, rate, standard_noise) 251 252 # Calculates the average teager value. 253 teager_value_of_std_wave = average_teager_value(standard_wave, amplitude) 254 255 return (teager_value_of_input / teager_value_of_std_wave) * standard_noise 256 257 258def error(f1, f2): 259 """Calculates the relative error between f1 and f2. 260 261 @param f1: Exact value. 262 @param f2: Test value. 263 264 @returns: Relative error between f1 and f2. 265 266 """ 267 return abs(float(f1) - float(f2)) / float(f1) 268 269 270def hilbert_analysis(signal, rate, block_size): 271 """Finds amplitude and frequency of each time of signal by Hilbert transform. 272 273 @param signal: The wave to analyze. 274 @param rate: Sampling rate. 275 @param block_size: The size of block to transform. 276 277 @returns: A tuple of list: (amplitude, frequency) composed of 278 amplitude and frequency of each time. 279 280 """ 281 # To apply Hilbert transform, the wave will be transformed 282 # segment by segment. For each segment, its size will be 283 # block_size and we will only take middle part of it. 284 # Thus, each segment looks like: |-----|=====|=====|-----|. 285 # "=...=" part will be taken while "-...-" part will be ignored. 286 # 287 # The whole size of taken part will be half of block_size 288 # which will be hilbert_block. 289 # The size of each ignored part will be half of hilbert_block 290 # which will be half_hilbert_block. 291 hilbert_block = block_size // 2 292 half_hilbert_block = hilbert_block // 2 293 # As mentioned above, for each block, we will only take middle 294 # part of it. Thus, the whole transformation will be completed as: 295 # |=====|=====|-----| |-----|=====|=====|-----| 296 # |-----|=====|=====|-----| |-----|=====|=====| 297 # |-----|=====|=====|-----| 298 # Specially, beginning and ending part may not have ignored part. 299 length = len(signal) 300 result = [] 301 for left_border in range(0, length, hilbert_block): 302 right_border = min(length, left_border + hilbert_block) 303 temp_left_border = max(0, left_border - half_hilbert_block) 304 temp_right_border = min(length, right_border + half_hilbert_block) 305 temp = hilbert(signal[temp_left_border:temp_right_border]) 306 for index in range(left_border, right_border): 307 result.append(temp[index - temp_left_border]) 308 result = numpy.asarray(result) 309 amplitude = numpy.abs(result) 310 phase = numpy.unwrap(numpy.angle(result)) 311 frequency = numpy.diff(phase) / (2.0 * numpy.pi) * rate 312 #frequency.append(frequency[len(frequency)-1]) 313 frequecny = numpy.append(frequency, frequency[len(frequency) - 1]) 314 return (amplitude, frequency) 315 316 317def find_block_average_value(arr, side_block_size, block_size): 318 """For each index, finds average value of its block, left block, right block. 319 320 It will find average value for each index in the range. 321 322 For each index, the range of its block is 323 [max(0, index - block_size / 2), min(length - 1, index + block_size / 2)] 324 For each index, the range of its left block is 325 [max(0, index - size_block_size), index] 326 For each index, the range of its right block is 327 [index, min(length - 1, index + side_block_size)] 328 329 @param arr: The array to be computed. 330 @param side_block_size: the size of the left_block and right_block. 331 @param block_size: the size of the block. 332 333 @returns: A tuple of lists: (left_block_average_array, 334 right_block_average_array, 335 block_average_array) 336 """ 337 length = len(arr) 338 left_border, right_border = 0, 1 339 left_block_sum = arr[0] 340 right_block_sum = arr[0] 341 left_average_array = numpy.zeros(length) 342 right_average_array = numpy.zeros(length) 343 block_average_array = numpy.zeros(length) 344 for index in range(0, length): 345 while left_border < index - side_block_size: 346 left_block_sum -= arr[left_border] 347 left_border += 1 348 while right_border < min(length, index + side_block_size): 349 right_block_sum += arr[right_border] 350 right_border += 1 351 352 left_average_value = float(left_block_sum) / (index - left_border + 1) 353 right_average_value = float(right_block_sum) / (right_border - index) 354 left_average_array[index] = left_average_value 355 right_average_array[index] = right_average_value 356 357 if index + 1 < length: 358 left_block_sum += arr[index + 1] 359 right_block_sum -= arr[index] 360 left_border, right_border = 0, 1 361 block_sum = 0 362 for index in range(0, length): 363 while left_border < index - block_size // 2: 364 block_sum -= arr[left_border] 365 left_border += 1 366 while right_border < min(length, index + block_size // 2): 367 block_sum += arr[right_border] 368 right_border += 1 369 370 average_value = float(block_sum) / (right_border - left_border) 371 block_average_array[index] = average_value 372 return (left_average_array, right_average_array, block_average_array) 373 374 375def find_start_end_index(dominant_frequency, 376 block_frequency_delta, 377 block_size, 378 frequency_error_threshold): 379 """Finds start and end index of sine wave. 380 381 For each block with size of block_size, we check that whether its frequency 382 is close enough to the dominant_frequency. If yes, we will consider this 383 block to be within the sine wave. 384 Then, it will return the start and end index of sine wave indicating that 385 sine wave is between [start_index, end_index) 386 It's okay if the whole signal only contains sine wave. 387 388 @param dominant_frequency: Dominant frequency of signal. 389 @param block_frequency_delta: Average absolute difference between dominant 390 frequency and frequency of each block. For 391 each index, its block is 392 [max(0, index - block_size / 2), 393 min(length - 1, index + block_size / 2)] 394 @param block_size: Block size in samples. 395 396 @returns: A tuple composed of (start_index, end_index) 397 398 """ 399 length = len(block_frequency_delta) 400 401 # Finds the start/end time index of playing based on dominant frequency 402 start_index, end_index = length - 1 , 0 403 for index in range(0, length): 404 left_border = max(0, index - block_size // 2) 405 right_border = min(length - 1, index + block_size // 2) 406 frequency_error = block_frequency_delta[index] / dominant_frequency 407 if frequency_error < frequency_error_threshold: 408 start_index = min(start_index, left_border) 409 end_index = max(end_index, right_border + 1) 410 return (start_index, end_index) 411 412 413def noise_detection(start_index, end_index, 414 block_amplitude, average_amplitude, 415 rate, noise_amplitude_threshold): 416 """Detects noise before/after sine wave. 417 418 If average amplitude of some sample's block before start of wave or after 419 end of wave is more than average_amplitude times noise_amplitude_threshold, 420 it will be considered as a noise. 421 422 @param start_index: Start index of sine wave. 423 @param end_index: End index of sine wave. 424 @param block_amplitude: An array for average amplitude of each block, where 425 amplitude is computed from Hilbert transform. 426 @param average_amplitude: Average amplitude of sine wave. 427 @param rate: Sampling rate 428 @param noise_amplitude_threshold: If the average amplitude of a block is 429 higher than average amplitude of the wave times 430 noise_amplitude_threshold, it will be considered as 431 noise before/after playback. 432 433 @returns: A tuple of lists indicating the time that noise happens: 434 (noise_before_playing, noise_after_playing). 435 436 """ 437 length = len(block_amplitude) 438 amplitude_threshold = average_amplitude * noise_amplitude_threshold 439 same_event_samples = rate * DEFAULT_SAME_EVENT_SECS 440 441 # Detects noise before playing. 442 noise_time_point = [] 443 last_noise_end_time_point = [] 444 previous_noise_index = None 445 times = 0 446 for index in range(0, length): 447 # Ignore noise too close to the beginning or the end of sine wave. 448 # Check the docstring of NEAR_SINE_START_OR_END_SECS. 449 if ((start_index - rate * NEAR_SINE_START_OR_END_SECS) <= index and 450 (index < end_index + rate * NEAR_SINE_START_OR_END_SECS)): 451 continue 452 453 # Ignore noise too close to the beginning or the end of original data. 454 # Check the docstring of NEAR_DATA_START_OR_END_SECS. 455 if (float(index) / rate <= 456 NEAR_DATA_START_OR_END_SECS + APPEND_ZEROS_SECS): 457 continue 458 if (float(length - index) / rate <= 459 NEAR_DATA_START_OR_END_SECS + APPEND_ZEROS_SECS): 460 continue 461 if block_amplitude[index] > amplitude_threshold: 462 same_event = False 463 if previous_noise_index: 464 same_event = (index - previous_noise_index) < same_event_samples 465 if not same_event: 466 index_start_sec = float(index) / rate - APPEND_ZEROS_SECS 467 index_end_sec = float(index + 1) / rate - APPEND_ZEROS_SECS 468 noise_time_point.append(index_start_sec) 469 last_noise_end_time_point.append(index_end_sec) 470 times += 1 471 index_end_sec = float(index + 1) / rate - APPEND_ZEROS_SECS 472 last_noise_end_time_point[times - 1] = index_end_sec 473 previous_noise_index = index 474 475 noise_before_playing, noise_after_playing = [], [] 476 for i in range(times): 477 duration = last_noise_end_time_point[i] - noise_time_point[i] 478 if noise_time_point[i] < float(start_index) / rate - APPEND_ZEROS_SECS: 479 noise_before_playing.append((noise_time_point[i], duration)) 480 else: 481 noise_after_playing.append((noise_time_point[i], duration)) 482 483 return (noise_before_playing, noise_after_playing) 484 485 486def delay_detection(start_index, end_index, 487 block_amplitude, average_amplitude, 488 dominant_frequency, 489 rate, 490 left_block_amplitude, 491 right_block_amplitude, 492 block_frequency_delta, 493 delay_amplitude_threshold, 494 frequency_error_threshold): 495 """Detects delay during playing. 496 497 For each sample, we will check whether the average amplitude of its block 498 is less than average amplitude of its left block and its right block times 499 delay_amplitude_threshold. Also, we will check whether the frequency of 500 its block is far from the dominant frequency. 501 If at least one constraint fulfilled, it will be considered as a delay. 502 503 @param start_index: Start index of sine wave. 504 @param end_index: End index of sine wave. 505 @param block_amplitude: An array for average amplitude of each block, where 506 amplitude is computed from Hilbert transform. 507 @param average_amplitude: Average amplitude of sine wave. 508 @param dominant_frequency: Dominant frequency of signal. 509 @param rate: Sampling rate 510 @param left_block_amplitude: Average amplitude of left block of each index. 511 Ref to find_block_average_value function. 512 @param right_block_amplitude: Average amplitude of right block of each index. 513 Ref to find_block_average_value function. 514 @param block_frequency_delta: Average absolute difference frequency to 515 dominant frequency of block of each index. 516 Ref to find_block_average_value function. 517 @param delay_amplitude_threshold: If the average amplitude of a block is 518 lower than average amplitude of the wave times 519 delay_amplitude_threshold, it will be considered 520 as delay. 521 @param frequency_error_threshold: Ref to DEFAULT_FREQUENCY_ERROR 522 523 @returns: List of delay occurrence: 524 [(time_1, duration_1), (time_2, duration_2), ...], 525 where time and duration are in seconds. 526 527 """ 528 delay_time_points = [] 529 last_delay_end_time_points = [] 530 previous_delay_index = None 531 times = 0 532 same_event_samples = rate * DEFAULT_SAME_EVENT_SECS 533 start_time = float(start_index) / rate - APPEND_ZEROS_SECS 534 end_time = float(end_index) / rate - APPEND_ZEROS_SECS 535 for index in range(start_index, end_index): 536 if block_amplitude[index] > average_amplitude * delay_amplitude_threshold: 537 continue 538 now_time = float(index) / rate - APPEND_ZEROS_SECS 539 if abs(now_time - start_time) < NEAR_START_OR_END_SECS: 540 continue 541 if abs(now_time - end_time) < NEAR_START_OR_END_SECS: 542 continue 543 # If amplitude less than its left/right side and small enough, 544 # it will be considered as a delay. 545 amp_threshold = average_amplitude * delay_amplitude_threshold 546 left_threshold = delay_amplitude_threshold * left_block_amplitude[index] 547 amp_threshold = min(amp_threshold, left_threshold) 548 right_threshold = delay_amplitude_threshold * right_block_amplitude[index] 549 amp_threshold = min(amp_threshold, right_threshold) 550 551 frequency_error = block_frequency_delta[index] / dominant_frequency 552 553 amplitude_too_small = block_amplitude[index] < amp_threshold 554 frequency_not_match = frequency_error > frequency_error_threshold 555 556 if amplitude_too_small or frequency_not_match: 557 same_event = False 558 if previous_delay_index: 559 same_event = (index - previous_delay_index) < same_event_samples 560 if not same_event: 561 index_start_sec = float(index) / rate - APPEND_ZEROS_SECS 562 index_end_sec = float(index + 1) / rate - APPEND_ZEROS_SECS 563 delay_time_points.append(index_start_sec) 564 last_delay_end_time_points.append(index_end_sec) 565 times += 1 566 previous_delay_index = index 567 index_end_sec = float(index + 1) / rate - APPEND_ZEROS_SECS 568 last_delay_end_time_points[times - 1] = index_end_sec 569 570 delay_list = [] 571 for i in range(len(delay_time_points)): 572 duration = last_delay_end_time_points[i] - delay_time_points[i] 573 delay_list.append( (delay_time_points[i], duration) ) 574 return delay_list 575 576 577def burst_detection(start_index, end_index, 578 block_amplitude, average_amplitude, 579 dominant_frequency, 580 rate, 581 left_block_amplitude, 582 right_block_amplitude, 583 block_frequency_delta, 584 burst_amplitude_threshold, 585 frequency_error_threshold): 586 """Detects burst during playing. 587 588 For each sample, we will check whether the average amplitude of its block is 589 more than average amplitude of its left block and its right block times 590 burst_amplitude_threshold. Also, we will check whether the frequency of 591 its block is not compatible to the dominant frequency. 592 If at least one constraint fulfilled, it will be considered as a burst. 593 594 @param start_index: Start index of sine wave. 595 @param end_index: End index of sine wave. 596 @param block_amplitude: An array for average amplitude of each block, where 597 amplitude is computed from Hilbert transform. 598 @param average_amplitude: Average amplitude of sine wave. 599 @param dominant_frequency: Dominant frequency of signal. 600 @param rate: Sampling rate 601 @param left_block_amplitude: Average amplitude of left block of each index. 602 Ref to find_block_average_value function. 603 @param right_block_amplitude: Average amplitude of right block of each index. 604 Ref to find_block_average_value function. 605 @param block_frequency_delta: Average absolute difference frequency to 606 dominant frequency of block of each index. 607 @param burst_amplitude_threshold: If the amplitude is higher than average 608 amplitude of its left block and its right block 609 times burst_amplitude_threshold. It will be 610 considered as a burst. 611 @param frequency_error_threshold: Ref to DEFAULT_FREQUENCY_ERROR 612 613 @returns: List of burst occurence: [time_1, time_2, ...], 614 where time is in seconds. 615 616 """ 617 burst_time_points = [] 618 previous_burst_index = None 619 same_event_samples = rate * DEFAULT_SAME_EVENT_SECS 620 for index in range(start_index, end_index): 621 # If amplitude higher than its left/right side and large enough, 622 # it will be considered as a burst. 623 if block_amplitude[index] <= average_amplitude * DEFAULT_BURST_TOO_SMALL: 624 continue 625 if abs(index - start_index) < rate * NEAR_START_OR_END_SECS: 626 continue 627 if abs(index - end_index) < rate * NEAR_START_OR_END_SECS: 628 continue 629 amp_threshold = average_amplitude * DEFAULT_BURST_TOO_SMALL 630 left_threshold = burst_amplitude_threshold * left_block_amplitude[index] 631 amp_threshold = max(amp_threshold, left_threshold) 632 right_threshold = burst_amplitude_threshold * right_block_amplitude[index] 633 amp_threshold = max(amp_threshold, right_threshold) 634 635 frequency_error = block_frequency_delta[index] / dominant_frequency 636 637 amplitude_too_large = block_amplitude[index] > amp_threshold 638 frequency_not_match = frequency_error > frequency_error_threshold 639 640 if amplitude_too_large or frequency_not_match: 641 same_event = False 642 if previous_burst_index: 643 same_event = index - previous_burst_index < same_event_samples 644 if not same_event: 645 burst_time_points.append(float(index) / rate - APPEND_ZEROS_SECS) 646 previous_burst_index = index 647 648 return burst_time_points 649 650 651def changing_volume_detection(start_index, end_index, 652 average_amplitude, 653 rate, 654 left_block_amplitude, 655 right_block_amplitude, 656 volume_changing_amplitude_threshold): 657 """Finds volume changing during playback. 658 659 For each index, we will compare average amplitude of its left block and its 660 right block. If average amplitude of right block is more than average 661 amplitude of left block times (1 + DEFAULT_VOLUME_CHANGE_AMPLITUDE), it will 662 be considered as an increasing volume. If the one of right block is less 663 than that of left block times (1 - DEFAULT_VOLUME_CHANGE_AMPLITUDE), it will 664 be considered as a decreasing volume. 665 666 @param start_index: Start index of sine wave. 667 @param end_index: End index of sine wave. 668 @param average_amplitude: Average amplitude of sine wave. 669 @param rate: Sampling rate 670 @param left_block_amplitude: Average amplitude of left block of each index. 671 Ref to find_block_average_value function. 672 @param right_block_amplitude: Average amplitude of right block of each index. 673 Ref to find_block_average_value function. 674 @param volume_changing_amplitude_threshold: If the average amplitude of right 675 block is higher or lower than 676 that of left one times this 677 value, it will be considered as 678 a volume change. 679 Also refer to 680 DEFAULT_VOLUME_CHANGE_AMPLITUDE 681 682 @returns: List of volume changing composed of 1 for increasing and 683 -1 for decreasing. 684 685 """ 686 length = len(left_block_amplitude) 687 688 # Detects rising and/or falling volume. 689 previous_rising_index, previous_falling_index = None, None 690 changing_time = [] 691 changing_events = [] 692 amplitude_threshold = average_amplitude * DEFAULT_VOLUME_CHANGE_TOO_SMALL 693 same_event_samples = rate * DEFAULT_SAME_EVENT_SECS 694 for index in range(start_index, end_index): 695 # Skips if amplitude is too small. 696 if left_block_amplitude[index] < amplitude_threshold: 697 continue 698 if right_block_amplitude[index] < amplitude_threshold: 699 continue 700 # Skips if changing is from start or end time 701 if float(abs(start_index - index)) / rate < NEAR_START_OR_END_SECS: 702 continue 703 if float(abs(end_index - index)) / rate < NEAR_START_OR_END_SECS: 704 continue 705 706 delta_margin = volume_changing_amplitude_threshold 707 if left_block_amplitude[index] > 0: 708 delta_margin *= left_block_amplitude[index] 709 710 increasing_threshold = left_block_amplitude[index] + delta_margin 711 decreasing_threshold = left_block_amplitude[index] - delta_margin 712 713 if right_block_amplitude[index] > increasing_threshold: 714 same_event = False 715 if previous_rising_index: 716 same_event = index - previous_rising_index < same_event_samples 717 if not same_event: 718 changing_time.append(float(index) / rate - APPEND_ZEROS_SECS) 719 changing_events.append(+1) 720 previous_rising_index = index 721 if right_block_amplitude[index] < decreasing_threshold: 722 same_event = False 723 if previous_falling_index: 724 same_event = index - previous_falling_index < same_event_samples 725 if not same_event: 726 changing_time.append(float(index) / rate - APPEND_ZEROS_SECS) 727 changing_events.append(-1) 728 previous_falling_index = index 729 730 # Combines consecutive increasing/decreasing event. 731 combined_changing_events, prev = [], 0 732 for i in range(len(changing_events)): 733 if changing_events[i] == prev: 734 continue 735 combined_changing_events.append((changing_time[i], changing_events[i])) 736 prev = changing_events[i] 737 return combined_changing_events 738 739 740def quality_measurement( 741 signal, rate, 742 dominant_frequency=None, 743 block_size_secs=DEFAULT_BLOCK_SIZE_SECS, 744 frequency_error_threshold=DEFAULT_FREQUENCY_ERROR, 745 delay_amplitude_threshold=DEFAULT_DELAY_AMPLITUDE_THRESHOLD, 746 noise_amplitude_threshold=DEFAULT_NOISE_AMPLITUDE_THRESHOLD, 747 burst_amplitude_threshold=DEFAULT_BURST_AMPLITUDE_THRESHOLD, 748 volume_changing_amplitude_threshold=DEFAULT_VOLUME_CHANGE_AMPLITUDE): 749 """Detects several artifacts and estimates the noise level. 750 751 This method detects artifact before playing, after playing, and delay 752 during playing. Also, it estimates the noise level of the signal. 753 To avoid the influence of noise, it calculates amplitude and frequency 754 block by block. 755 756 @param signal: A list of numbers for one-channel PCM data. The data should 757 be normalized to [-1,1]. 758 @param rate: Sampling rate 759 @param dominant_frequency: Dominant frequency of signal. Set None to 760 recalculate the frequency in this function. 761 @param block_size_secs: Block size in seconds. The measurement will be done 762 block-by-block using average amplitude and frequency 763 in each block to avoid noise. 764 @param frequency_error_threshold: Ref to DEFAULT_FREQUENCY_ERROR. 765 @param delay_amplitude_threshold: If the average amplitude of a block is 766 lower than average amplitude of the wave 767 times delay_amplitude_threshold, it will 768 be considered as delay. 769 Also refer to delay_detection and 770 DEFAULT_DELAY_AMPLITUDE_THRESHOLD. 771 @param noise_amplitude_threshold: If the average amplitude of a block is 772 higher than average amplitude of the wave 773 times noise_amplitude_threshold, it will 774 be considered as noise before/after 775 playback. 776 Also refer to noise_detection and 777 DEFAULT_NOISE_AMPLITUDE_THRESHOLD. 778 @param burst_amplitude_threshold: If the average amplitude of a block is 779 higher than average amplitude of its left 780 block and its right block times 781 burst_amplitude_threshold. It will be 782 considered as a burst. 783 Also refer to burst_detection and 784 DEFAULT_BURST_AMPLITUDE_THRESHOLD. 785 @param volume_changing_amplitude_threshold: If the average amplitude of right 786 block is higher or lower than 787 that of left one times this 788 value, it will be considered as 789 a volume change. 790 Also refer to 791 changing_volume_detection and 792 DEFAULT_VOLUME_CHANGE_AMPLITUDE 793 794 @returns: A dictoinary of detection/estimation: 795 {'artifacts': 796 {'noise_before_playback': 797 [(time_1, duration_1), (time_2, duration_2), ...], 798 'noise_after_playback': 799 [(time_1, duration_1), (time_2, duration_2), ...], 800 'delay_during_playback': 801 [(time_1, duration_1), (time_2, duration_2), ...], 802 'burst_during_playback': 803 [time_1, time_2, ...] 804 }, 805 'volume_changes': 806 [(time_1, flag_1), (time_2, flag_2), ...], 807 'equivalent_noise_level': level 808 } 809 where durations and time points are in seconds. And, 810 equivalence_noise_level is the quotient of noise and 811 wave which refers to DEFAULT_STANDARD_NOISE. 812 volume_changes is a list of tuples containing time 813 stamps and decreasing/increasing flags for volume 814 change events. 815 816 """ 817 # Calculates the block size, from seconds to samples. 818 block_size = int(block_size_secs * rate) 819 820 signal = numpy.concatenate((numpy.zeros(int(rate * APPEND_ZEROS_SECS)), 821 signal, 822 numpy.zeros(int(rate * APPEND_ZEROS_SECS)))) 823 signal = numpy.array(signal, dtype=float) 824 length = len(signal) 825 826 # Calculates the amplitude and frequency. 827 amplitude, frequency = hilbert_analysis(signal, rate, block_size) 828 829 # Finds the dominant frequency. 830 if not dominant_frequency: 831 dominant_frequency = audio_analysis.spectral_analysis(signal, rate)[0][0] 832 833 # Finds the array which contains absolute difference between dominant 834 # frequency and frequency at each time point. 835 frequency_delta = abs(frequency - dominant_frequency) 836 837 # Computes average amplitude of each type of block 838 res = find_block_average_value(amplitude, block_size * 2, block_size) 839 left_block_amplitude, right_block_amplitude, block_amplitude = res 840 841 # Computes average absolute difference of frequency and dominant frequency 842 # of the block of each index 843 _, _, block_frequency_delta = find_block_average_value(frequency_delta, 844 block_size * 2, 845 block_size) 846 847 # Finds start and end index of sine wave. 848 start_index, end_index = find_start_end_index(dominant_frequency, 849 block_frequency_delta, 850 block_size, 851 frequency_error_threshold) 852 853 if start_index > end_index: 854 raise SineWaveNotFound('No sine wave found in signal') 855 856 logging.debug('Found sine wave: start: %s, end: %s', 857 float(start_index) / rate - APPEND_ZEROS_SECS, 858 float(end_index) / rate - APPEND_ZEROS_SECS) 859 860 sum_of_amplitude = float(sum(amplitude[start_index:end_index])) 861 # Finds average amplitude of sine wave. 862 average_amplitude = sum_of_amplitude / (end_index - start_index) 863 864 # Finds noise before and/or after playback. 865 noise_before_playing, noise_after_playing = noise_detection( 866 start_index, end_index, 867 block_amplitude, average_amplitude, 868 rate, 869 noise_amplitude_threshold) 870 871 # Finds delay during playback. 872 delays = delay_detection(start_index, end_index, 873 block_amplitude, average_amplitude, 874 dominant_frequency, 875 rate, 876 left_block_amplitude, 877 right_block_amplitude, 878 block_frequency_delta, 879 delay_amplitude_threshold, 880 frequency_error_threshold) 881 882 # Finds burst during playback. 883 burst_time_points = burst_detection(start_index, end_index, 884 block_amplitude, average_amplitude, 885 dominant_frequency, 886 rate, 887 left_block_amplitude, 888 right_block_amplitude, 889 block_frequency_delta, 890 burst_amplitude_threshold, 891 frequency_error_threshold) 892 893 # Finds volume changing during playback. 894 volume_changes = changing_volume_detection( 895 start_index, end_index, 896 average_amplitude, 897 rate, 898 left_block_amplitude, 899 right_block_amplitude, 900 volume_changing_amplitude_threshold) 901 902 # Calculates the average teager value. 903 teager_value = average_teager_value(signal[start_index:end_index], 904 average_amplitude) 905 906 # Finds out the noise level. 907 noise = noise_level(average_amplitude, dominant_frequency, 908 rate, 909 teager_value) 910 911 return {'artifacts': 912 {'noise_before_playback': noise_before_playing, 913 'noise_after_playback': noise_after_playing, 914 'delay_during_playback': delays, 915 'burst_during_playback': burst_time_points 916 }, 917 'volume_changes': volume_changes, 918 'equivalent_noise_level': noise 919 } 920