1 /* 2 * Copyright (C) 2017 The Android Open Source Project 3 * 4 * Licensed under the Apache License, Version 2.0 (the "License"); 5 * you may not use this file except in compliance with the License. 6 * You may obtain a copy of the License at 7 * 8 * http://www.apache.org/licenses/LICENSE-2.0 9 * 10 * Unless required by applicable law or agreed to in writing, software 11 * distributed under the License is distributed on an "AS IS" BASIS, 12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 13 * See the License for the specific language governing permissions and 14 * limitations under the License. 15 */ 16 17 #ifndef ANALYZER_GLITCH_ANALYZER_H 18 #define ANALYZER_GLITCH_ANALYZER_H 19 20 #include <algorithm> 21 #include <cctype> 22 #include <iomanip> 23 #include <iostream> 24 25 #include "LatencyAnalyzer.h" 26 #include "PseudoRandom.h" 27 28 /** 29 * Output a steady sine wave and analyze the return signal. 30 * 31 * Use a cosine transform to measure the predicted magnitude and relative phase of the 32 * looped back sine wave. Then generate a predicted signal and compare with the actual signal. 33 */ 34 class GlitchAnalyzer : public LoopbackProcessor { 35 public: 36 getState()37 int32_t getState() const { 38 return mState; 39 } 40 getPeakAmplitude()41 double getPeakAmplitude() const { 42 return mPeakFollower.getLevel(); 43 } 44 getTolerance()45 double getTolerance() { 46 return mTolerance; 47 } 48 setTolerance(double tolerance)49 void setTolerance(double tolerance) { 50 mTolerance = tolerance; 51 mScaledTolerance = mMagnitude * mTolerance; 52 } 53 setMagnitude(double magnitude)54 void setMagnitude(double magnitude) { 55 mMagnitude = magnitude; 56 mScaledTolerance = mMagnitude * mTolerance; 57 } 58 getGlitchCount()59 int32_t getGlitchCount() const { 60 return mGlitchCount; 61 } 62 getStateFrameCount(int state)63 int32_t getStateFrameCount(int state) const { 64 return mStateFrameCounters[state]; 65 } 66 getSignalToNoiseDB()67 double getSignalToNoiseDB() { 68 static const double threshold = 1.0e-14; 69 if (mMeanSquareSignal < threshold || mMeanSquareNoise < threshold) { 70 return 0.0; 71 } else { 72 double signalToNoise = mMeanSquareSignal / mMeanSquareNoise; // power ratio 73 double signalToNoiseDB = 10.0 * log(signalToNoise); 74 if (signalToNoiseDB < MIN_SNR_DB) { 75 ALOGD("ERROR - signal to noise ratio is too low! < %d dB. Adjust volume.", 76 MIN_SNR_DB); 77 setResult(ERROR_VOLUME_TOO_LOW); 78 } 79 return signalToNoiseDB; 80 } 81 } 82 analyze()83 std::string analyze() override { 84 std::stringstream report; 85 report << "GlitchAnalyzer ------------------\n"; 86 report << LOOPBACK_RESULT_TAG "peak.amplitude = " << std::setw(8) 87 << getPeakAmplitude() << "\n"; 88 report << LOOPBACK_RESULT_TAG "sine.magnitude = " << std::setw(8) 89 << mMagnitude << "\n"; 90 report << LOOPBACK_RESULT_TAG "rms.noise = " << std::setw(8) 91 << mMeanSquareNoise << "\n"; 92 report << LOOPBACK_RESULT_TAG "signal.to.noise.db = " << std::setw(8) 93 << getSignalToNoiseDB() << "\n"; 94 report << LOOPBACK_RESULT_TAG "frames.accumulated = " << std::setw(8) 95 << mFramesAccumulated << "\n"; 96 report << LOOPBACK_RESULT_TAG "sine.period = " << std::setw(8) 97 << mSinePeriod << "\n"; 98 report << LOOPBACK_RESULT_TAG "test.state = " << std::setw(8) 99 << mState << "\n"; 100 report << LOOPBACK_RESULT_TAG "frame.count = " << std::setw(8) 101 << mFrameCounter << "\n"; 102 // Did we ever get a lock? 103 bool gotLock = (mState == STATE_LOCKED) || (mGlitchCount > 0); 104 if (!gotLock) { 105 report << "ERROR - failed to lock on reference sine tone.\n"; 106 setResult(ERROR_NO_LOCK); 107 } else { 108 // Only print if meaningful. 109 report << LOOPBACK_RESULT_TAG "glitch.count = " << std::setw(8) 110 << mGlitchCount << "\n"; 111 report << LOOPBACK_RESULT_TAG "max.glitch = " << std::setw(8) 112 << mMaxGlitchDelta << "\n"; 113 if (mGlitchCount > 0) { 114 report << "ERROR - number of glitches > 0\n"; 115 setResult(ERROR_GLITCHES); 116 } 117 } 118 return report.str(); 119 } 120 printStatus()121 void printStatus() override { 122 ALOGD("st = %d, #gl = %3d,", mState, mGlitchCount); 123 } 124 /** 125 * Calculate the magnitude of the component of the input signal 126 * that matches the analysis frequency. 127 * Also calculate the phase that we can use to create a 128 * signal that matches that component. 129 * The phase will be between -PI and +PI. 130 */ 131 double calculateMagnitude(double *phasePtr = nullptr) { 132 if (mFramesAccumulated == 0) { 133 return 0.0; 134 } 135 double sinMean = mSinAccumulator / mFramesAccumulated; 136 double cosMean = mCosAccumulator / mFramesAccumulated; 137 double magnitude = 2.0 * sqrt((sinMean * sinMean) + (cosMean * cosMean)); 138 if (phasePtr != nullptr) { 139 double phase = M_PI_2 - atan2(sinMean, cosMean); 140 *phasePtr = phase; 141 } 142 return magnitude; 143 } 144 145 /** 146 * @param frameData contains microphone data with sine signal feedback 147 * @param channelCount 148 */ processInputFrame(float * frameData,int)149 result_code processInputFrame(float *frameData, int /* channelCount */) override { 150 result_code result = RESULT_OK; 151 152 float sample = frameData[0]; 153 float peak = mPeakFollower.process(sample); 154 155 // Force a periodic glitch to test the detector! 156 if (mForceGlitchDuration > 0) { 157 if (mForceGlitchCounter == 0) { 158 ALOGE("%s: force a glitch!!", __func__); 159 mForceGlitchCounter = getSampleRate(); 160 } else if (mForceGlitchCounter <= mForceGlitchDuration) { 161 // Force an abrupt offset. 162 sample += (sample > 0.0) ? -0.5f : 0.5f; 163 } 164 --mForceGlitchCounter; 165 } 166 167 mStateFrameCounters[mState]++; // count how many frames we are in each state 168 169 switch (mState) { 170 case STATE_IDLE: 171 mDownCounter--; 172 if (mDownCounter <= 0) { 173 mState = STATE_IMMUNE; 174 mDownCounter = IMMUNE_FRAME_COUNT; 175 mInputPhase = 0.0; // prevent spike at start 176 mOutputPhase = 0.0; 177 } 178 break; 179 180 case STATE_IMMUNE: 181 mDownCounter--; 182 if (mDownCounter <= 0) { 183 mState = STATE_WAITING_FOR_SIGNAL; 184 } 185 break; 186 187 case STATE_WAITING_FOR_SIGNAL: 188 if (peak > mThreshold) { 189 mState = STATE_WAITING_FOR_LOCK; 190 //ALOGD("%5d: switch to STATE_WAITING_FOR_LOCK", mFrameCounter); 191 resetAccumulator(); 192 } 193 break; 194 195 case STATE_WAITING_FOR_LOCK: 196 mSinAccumulator += sample * sinf(mInputPhase); 197 mCosAccumulator += sample * cosf(mInputPhase); 198 mFramesAccumulated++; 199 // Must be a multiple of the period or the calculation will not be accurate. 200 if (mFramesAccumulated == mSinePeriod * PERIODS_NEEDED_FOR_LOCK) { 201 double phaseOffset = 0.0; 202 setMagnitude(calculateMagnitude(&phaseOffset)); 203 // ALOGD("%s() mag = %f, offset = %f, prev = %f", 204 // __func__, mMagnitude, mPhaseOffset, mPreviousPhaseOffset); 205 if (mMagnitude > mThreshold) { 206 if (abs(phaseOffset) < kMaxPhaseError) { 207 mState = STATE_LOCKED; 208 // ALOGD("%5d: switch to STATE_LOCKED", mFrameCounter); 209 } 210 // Adjust mInputPhase to match measured phase 211 mInputPhase += phaseOffset; 212 } 213 resetAccumulator(); 214 } 215 incrementInputPhase(); 216 break; 217 218 case STATE_LOCKED: { 219 // Predict next sine value 220 double predicted = sinf(mInputPhase) * mMagnitude; 221 double diff = predicted - sample; 222 double absDiff = fabs(diff); 223 mMaxGlitchDelta = std::max(mMaxGlitchDelta, absDiff); 224 if (absDiff > mScaledTolerance) { 225 result = ERROR_GLITCHES; 226 onGlitchStart(); 227 // LOGI("diff glitch detected, absDiff = %g", absDiff); 228 } else { 229 mSumSquareSignal += predicted * predicted; 230 mSumSquareNoise += diff * diff; 231 // Track incoming signal and slowly adjust magnitude to account 232 // for drift in the DRC or AGC. 233 mSinAccumulator += sample * sinf(mInputPhase); 234 mCosAccumulator += sample * cosf(mInputPhase); 235 mFramesAccumulated++; 236 // Must be a multiple of the period or the calculation will not be accurate. 237 if (mFramesAccumulated == mSinePeriod) { 238 const double coefficient = 0.1; 239 double phaseOffset = 0.0; 240 double magnitude = calculateMagnitude(&phaseOffset); 241 // One pole averaging filter. 242 setMagnitude((mMagnitude * (1.0 - coefficient)) + (magnitude * coefficient)); 243 244 mMeanSquareNoise = mSumSquareNoise * mInverseSinePeriod; 245 mMeanSquareSignal = mSumSquareSignal * mInverseSinePeriod; 246 resetAccumulator(); 247 248 if (abs(phaseOffset) > kMaxPhaseError) { 249 result = ERROR_GLITCHES; 250 onGlitchStart(); 251 ALOGD("phase glitch detected, phaseOffset = %g", phaseOffset); 252 } else if (mMagnitude < mThreshold) { 253 result = ERROR_GLITCHES; 254 onGlitchStart(); 255 ALOGD("magnitude glitch detected, mMagnitude = %g", mMagnitude); 256 } 257 } 258 } 259 incrementInputPhase(); 260 } break; 261 262 case STATE_GLITCHING: { 263 // Predict next sine value 264 mGlitchLength++; 265 double predicted = sinf(mInputPhase) * mMagnitude; 266 double diff = predicted - sample; 267 double absDiff = fabs(diff); 268 mMaxGlitchDelta = std::max(mMaxGlitchDelta, absDiff); 269 if (absDiff < mScaledTolerance) { // close enough? 270 // If we get a full sine period of non-glitch samples in a row then consider the glitch over. 271 // We don't want to just consider a zero crossing the end of a glitch. 272 if (mNonGlitchCount++ > mSinePeriod) { 273 onGlitchEnd(); 274 } 275 } else { 276 mNonGlitchCount = 0; 277 if (mGlitchLength > (4 * mSinePeriod)) { 278 relock(); 279 } 280 } 281 incrementInputPhase(); 282 } break; 283 284 case NUM_STATES: // not a real state 285 break; 286 } 287 288 mFrameCounter++; 289 290 return result; 291 } 292 293 // advance and wrap phase incrementInputPhase()294 void incrementInputPhase() { 295 mInputPhase += mPhaseIncrement; 296 if (mInputPhase > M_PI) { 297 mInputPhase -= (2.0 * M_PI); 298 } 299 } 300 301 // advance and wrap phase incrementOutputPhase()302 void incrementOutputPhase() { 303 mOutputPhase += mPhaseIncrement; 304 if (mOutputPhase > M_PI) { 305 mOutputPhase -= (2.0 * M_PI); 306 } 307 } 308 309 /** 310 * @param frameData upon return, contains the reference sine wave 311 * @param channelCount 312 */ processOutputFrame(float * frameData,int channelCount)313 result_code processOutputFrame(float *frameData, int channelCount) override { 314 float output = 0.0f; 315 // Output sine wave so we can measure it. 316 if (mState != STATE_IDLE) { 317 float sinOut = sinf(mOutputPhase); 318 incrementOutputPhase(); 319 output = (sinOut * mOutputAmplitude) 320 + (mWhiteNoise.nextRandomDouble() * kNoiseAmplitude); 321 // ALOGD("sin(%f) = %f, %f\n", mOutputPhase, sinOut, mPhaseIncrement); 322 } 323 frameData[0] = output; 324 for (int i = 1; i < channelCount; i++) { 325 frameData[i] = 0.0f; 326 } 327 return RESULT_OK; 328 } 329 onGlitchStart()330 void onGlitchStart() { 331 mGlitchCount++; 332 // ALOGD("%5d: STARTED a glitch # %d", mFrameCounter, mGlitchCount); 333 mState = STATE_GLITCHING; 334 mGlitchLength = 1; 335 mNonGlitchCount = 0; 336 } 337 onGlitchEnd()338 void onGlitchEnd() { 339 // ALOGD("%5d: ENDED a glitch # %d, length = %d", mFrameCounter, mGlitchCount, mGlitchLength); 340 mState = STATE_LOCKED; 341 resetAccumulator(); 342 } 343 344 // reset the sine wave detector resetAccumulator()345 void resetAccumulator() { 346 mFramesAccumulated = 0; 347 mSinAccumulator = 0.0; 348 mCosAccumulator = 0.0; 349 mSumSquareSignal = 0.0; 350 mSumSquareNoise = 0.0; 351 } 352 relock()353 void relock() { 354 // ALOGD("relock: %d because of a very long %d glitch", mFrameCounter, mGlitchLength); 355 mState = STATE_WAITING_FOR_LOCK; 356 resetAccumulator(); 357 } 358 reset()359 void reset() override { 360 LoopbackProcessor::reset(); 361 mState = STATE_IDLE; 362 mDownCounter = IDLE_FRAME_COUNT; 363 resetAccumulator(); 364 } 365 prepareToTest()366 void prepareToTest() override { 367 LoopbackProcessor::prepareToTest(); 368 mSinePeriod = getSampleRate() / kTargetGlitchFrequency; 369 mOutputPhase = 0.0f; 370 mInverseSinePeriod = 1.0 / mSinePeriod; 371 mPhaseIncrement = 2.0 * M_PI * mInverseSinePeriod; 372 mGlitchCount = 0; 373 mMaxGlitchDelta = 0.0; 374 for (int i = 0; i < NUM_STATES; i++) { 375 mStateFrameCounters[i] = 0; 376 } 377 } 378 379 private: 380 381 // These must match the values in GlitchActivity.java 382 enum sine_state_t { 383 STATE_IDLE, // beginning 384 STATE_IMMUNE, // ignoring input, waiting fo HW to settle 385 STATE_WAITING_FOR_SIGNAL, // looking for a loud signal 386 STATE_WAITING_FOR_LOCK, // trying to lock onto the phase of the sine 387 STATE_LOCKED, // locked on the sine wave, looking for glitches 388 STATE_GLITCHING, // locked on the sine wave but glitching 389 NUM_STATES 390 }; 391 392 enum constants { 393 // Arbitrary durations, assuming 48000 Hz 394 IDLE_FRAME_COUNT = 48 * 100, 395 IMMUNE_FRAME_COUNT = 48 * 100, 396 PERIODS_NEEDED_FOR_LOCK = 8, 397 MIN_SNR_DB = 65 398 }; 399 400 static constexpr float kNoiseAmplitude = 0.00; // Used to experiment with warbling caused by DRC. 401 static constexpr int kTargetGlitchFrequency = 607; 402 static constexpr double kMaxPhaseError = M_PI * 0.05; 403 404 float mTolerance = 0.10; // scaled from 0.0 to 1.0 405 double mThreshold = 0.005; 406 int mSinePeriod = 1; // this will be set before use 407 double mInverseSinePeriod = 1.0; 408 409 int32_t mStateFrameCounters[NUM_STATES]; 410 411 double mPhaseIncrement = 0.0; 412 double mInputPhase = 0.0; 413 double mOutputPhase = 0.0; 414 double mMagnitude = 0.0; 415 int32_t mFramesAccumulated = 0; 416 double mSinAccumulator = 0.0; 417 double mCosAccumulator = 0.0; 418 double mMaxGlitchDelta = 0.0; 419 int32_t mGlitchCount = 0; 420 int32_t mNonGlitchCount = 0; 421 int32_t mGlitchLength = 0; 422 // This is used for processing every frame so we cache it here. 423 double mScaledTolerance = 0.0; 424 int mDownCounter = IDLE_FRAME_COUNT; 425 int32_t mFrameCounter = 0; 426 double mOutputAmplitude = 0.75; 427 428 int32_t mForceGlitchDuration = 0; // if > 0 then force a glitch for debugging 429 int32_t mForceGlitchCounter = 4 * 48000; // count down and trigger at zero 430 431 // measure background noise continuously as a deviation from the expected signal 432 double mSumSquareSignal = 0.0; 433 double mSumSquareNoise = 0.0; 434 double mMeanSquareSignal = 0.0; 435 double mMeanSquareNoise = 0.0; 436 437 PeakDetector mPeakFollower; 438 439 PseudoRandom mWhiteNoise; 440 441 sine_state_t mState = STATE_IDLE; 442 }; 443 444 445 #endif //ANALYZER_GLITCH_ANALYZER_H 446