drivers/staging/echo/echo.c

   1 /*
   2  * SpanDSP - a series of DSP components for telephony
   3  *
   4  * echo.c - A line echo canceller.  This code is being developed
   5  *          against and partially complies with G168.
   6  *
   7  * Written by Steve Underwood <steveu@coppice.org>
   8  *         and David Rowe <david_at_rowetel_dot_com>
   9  *
  10  * Copyright (C) 2001, 2003 Steve Underwood, 2007 David Rowe
  11  *
  12  * Based on a bit from here, a bit from there, eye of toad, ear of
  13  * bat, 15 years of failed attempts by David and a few fried brain
  14  * cells.
  15  *
  16  * All rights reserved.
  17  *
  18  * This program is free software; you can redistribute it and/or modify
  19  * it under the terms of the GNU General Public License version 2, as
  20  * published by the Free Software Foundation.
  21  *
  22  * This program is distributed in the hope that it will be useful,
  23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  25  * GNU General Public License for more details.
  26  *
  27  * You should have received a copy of the GNU General Public License
  28  * along with this program; if not, write to the Free Software
  29  * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  30  *
  31  * $Id: echo.c,v 1.20 2006/12/01 18:00:48 steveu Exp $
  32  */
  33
  34 /*! \file */
  35
  36 /* Implementation Notes
  37    David Rowe
  38    April 2007
  39
  40    This code started life as Steve's NLMS algorithm with a tap
  41    rotation algorithm to handle divergence during double talk.  I
  42    added a Geigel Double Talk Detector (DTD) [2] and performed some
  43    G168 tests.  However I had trouble meeting the G168 requirements,
  44    especially for double talk - there were always cases where my DTD
  45    failed, for example where near end speech was under the 6dB
  46    threshold required for declaring double talk.
  47
  48    So I tried a two path algorithm [1], which has so far given better
  49    results.  The original tap rotation/Geigel algorithm is available
  50    in SVN http://svn.rowetel.com/software/oslec/tags/before_16bit.
  51    It's probably possible to make it work if some one wants to put some
  52    serious work into it.
  53
  54    At present no special treatment is provided for tones, which
  55    generally cause NLMS algorithms to diverge.  Initial runs of a
  56    subset of the G168 tests for tones (e.g ./echo_test 6) show the
  57    current algorithm is passing OK, which is kind of surprising.  The
  58    full set of tests needs to be performed to confirm this result.
  59
  60    One other interesting change is that I have managed to get the NLMS
  61    code to work with 16 bit coefficients, rather than the original 32
  62    bit coefficents.  This reduces the MIPs and storage required.
  63    I evaulated the 16 bit port using g168_tests.sh and listening tests
  64    on 4 real-world samples.
  65
  66    I also attempted the implementation of a block based NLMS update
  67    [2] but although this passes g168_tests.sh it didn't converge well
  68    on the real-world samples.  I have no idea why, perhaps a scaling
  69    problem.  The block based code is also available in SVN
  70    http://svn.rowetel.com/software/oslec/tags/before_16bit.  If this
  71    code can be debugged, it will lead to further reduction in MIPS, as
  72    the block update code maps nicely onto DSP instruction sets (it's a
  73    dot product) compared to the current sample-by-sample update.
  74
  75    Steve also has some nice notes on echo cancellers in echo.h
  76
  77
  78    References:
  79
  80    [1] Ochiai, Areseki, and Ogihara, "Echo Canceller with Two Echo
  81        Path Models", IEEE Transactions on communications, COM-25,
  82        No. 6, June
  83        1977.
  84        http://www.rowetel.com/images/echo/dual_path_paper.pdf
  85
  86    [2] The classic, very useful paper that tells you how to
  87        actually build a real world echo canceller:
  88          Messerschmitt, Hedberg, Cole, Haoui, Winship, "Digital Voice
  89          Echo Canceller with a TMS320020,
  90          http://www.rowetel.com/images/echo/spra129.pdf
  91
  92    [3] I have written a series of blog posts on this work, here is
  93        Part 1: http://www.rowetel.com/blog/?p=18
  94
  95    [4] The source code http://svn.rowetel.com/software/oslec/
  96
  97    [5] A nice reference on LMS filters:
  98          http://en.wikipedia.org/wiki/Least_mean_squares_filter
  99
 100    Credits:
 101
 102    Thanks to Steve Underwood, Jean-Marc Valin, and Ramakrishnan
 103    Muthukrishnan for their suggestions and email discussions.  Thanks
 104    also to those people who collected echo samples for me such as
 105    Mark, Pawel, and Pavel.
 106 */
 107
 108 #include <linux/kernel.h>       /* We're doing kernel work */
 109 #include <linux/module.h>
 110 #include <linux/kernel.h>
 111 #include <linux/slab.h>
 112 #define malloc(a) kmalloc((a), GFP_KERNEL)
 113 #define free(a) kfree(a)
 114
 115 #include "bit_operations.h"
 116 #include "echo.h"
 117
 118 #define MIN_TX_POWER_FOR_ADAPTION   64
 119 #define MIN_RX_POWER_FOR_ADAPTION   64
 120 #define DTD_HANGOVER               600     /* 600 samples, or 75ms     */
 121 #define DC_LOG2BETA                  3     /* log2() of DC filter Beta */
 122
 123 /*-----------------------------------------------------------------------*\
 124                                FUNCTIONS
 125 \*-----------------------------------------------------------------------*/
 126
 127 /* adapting coeffs using the traditional stochastic descent (N)LMS algorithm */
 128
 129
 130 #ifdef __bfin__
 131 static void __inline__ lms_adapt_bg(struct oslec_state *ec, int clean, int shift)
 132 {
 133     int i, j;
 134     int offset1;
 135     int offset2;
 136     int factor;
 137     int exp;
 138     int16_t *phist;
 139     int n;
 140
 141     if (shift > 0)
 142         factor = clean << shift;
 143     else
 144         factor = clean >> -shift;
 145
 146     /* Update the FIR taps */
 147
 148     offset2 = ec->curr_pos;
 149     offset1 = ec->taps - offset2;
 150     phist = &ec->fir_state_bg.history[offset2];
 151
 152     /* st: and en: help us locate the assembler in echo.s */
 153
 154     //asm("st:");
 155     n = ec->taps;
 156     for (i = 0, j = offset2;  i < n;  i++, j++)
 157     {
 158        exp = *phist++ * factor;
 159        ec->fir_taps16[1][i] += (int16_t) ((exp+(1<<14)) >> 15);
 160     }
 161     //asm("en:");
 162
 163     /* Note the asm for the inner loop above generated by Blackfin gcc
 164        4.1.1 is pretty good (note even parallel instructions used):
 165
 166         R0 = W [P0++] (X);
 167         R0 *= R2;
 168         R0 = R0 + R3 (NS) ||
 169         R1 = W [P1] (X) ||
 170         nop;
 171         R0 >>>= 15;
 172         R0 = R0 + R1;
 173         W [P1++] = R0;
 174
 175         A block based update algorithm would be much faster but the
 176         above can't be improved on much.  Every instruction saved in
 177         the loop above is 2 MIPs/ch!  The for loop above is where the
 178         Blackfin spends most of it's time - about 17 MIPs/ch measured
 179         with speedtest.c with 256 taps (32ms).  Write-back and
 180         Write-through cache gave about the same performance.
 181     */
 182 }
 183
 184 /*
 185    IDEAS for further optimisation of lms_adapt_bg():
 186
 187    1/ The rounding is quite costly.  Could we keep as 32 bit coeffs
 188    then make filter pluck the MS 16-bits of the coeffs when filtering?
 189    However this would lower potential optimisation of filter, as I
 190    think the dual-MAC architecture requires packed 16 bit coeffs.
 191
 192    2/ Block based update would be more efficient, as per comments above,
 193    could use dual MAC architecture.
 194
 195    3/ Look for same sample Blackfin LMS code, see if we can get dual-MAC
 196    packing.
 197
 198    4/ Execute the whole e/c in a block of say 20ms rather than sample
 199    by sample.  Processing a few samples every ms is inefficient.
 200 */
 201
 202 #else
 203 static __inline__ void lms_adapt_bg(struct oslec_state *ec, int clean, int shift)
 204 {
 205     int i;
 206
 207     int offset1;
 208     int offset2;
 209     int factor;
 210     int exp;
 211
 212     if (shift > 0)
 213         factor = clean << shift;
 214     else
 215         factor = clean >> -shift;
 216
 217     /* Update the FIR taps */
 218
 219     offset2 = ec->curr_pos;
 220     offset1 = ec->taps - offset2;
 221
 222     for (i = ec->taps - 1;  i >= offset1;  i--)
 223     {
 224        exp = (ec->fir_state_bg.history[i - offset1]*factor);
 225        ec->fir_taps16[1][i] += (int16_t) ((exp+(1<<14)) >> 15);
 226     }
 227     for (  ;  i >= 0;  i--)
 228     {
 229        exp = (ec->fir_state_bg.history[i + offset2]*factor);
 230        ec->fir_taps16[1][i] += (int16_t) ((exp+(1<<14)) >> 15);
 231     }
 232 }
 233 #endif
 234
 235 /*- End of function --------------------------------------------------------*/
 236
 237 struct oslec_state *oslec_create(int len, int adaption_mode)
 238 {
 239     struct oslec_state *ec;
 240     int i;
 241     int j;
 242
 243     ec = kmalloc(sizeof(*ec), GFP_KERNEL);
 244     if (ec == NULL)
 245         return  NULL;
 246     memset(ec, 0, sizeof(*ec));
 247
 248     ec->taps = len;
 249     ec->log2taps = top_bit(len);
 250     ec->curr_pos = ec->taps - 1;
 251
 252     for (i = 0;  i < 2;  i++)
 253     {
 254         if ((ec->fir_taps16[i] = (int16_t *) malloc((ec->taps)*sizeof(int16_t))) == NULL)
 255         {
 256             for (j = 0;  j < i;  j++)
 257                 kfree(ec->fir_taps16[j]);
 258             kfree(ec);
 259             return  NULL;
 260         }
 261         memset(ec->fir_taps16[i], 0, (ec->taps)*sizeof(int16_t));
 262     }
 263
 264     fir16_create(&ec->fir_state,
 265                  ec->fir_taps16[0],
 266                  ec->taps);
 267     fir16_create(&ec->fir_state_bg,
 268                  ec->fir_taps16[1],
 269                  ec->taps);
 270
 271     for(i=0; i<5; i++) {
 272       ec->xvtx[i] = ec->yvtx[i] = ec->xvrx[i] = ec->yvrx[i] = 0;
 273     }
 274
 275     ec->cng_level = 1000;
 276     oslec_adaption_mode(ec, adaption_mode);
 277
 278     ec->snapshot = (int16_t*)malloc(ec->taps*sizeof(int16_t));
 279     memset(ec->snapshot, 0, sizeof(int16_t)*ec->taps);
 280
 281     ec->cond_met = 0;
 282     ec->Pstates = 0;
 283     ec->Ltxacc = ec->Lrxacc = ec->Lcleanacc = ec->Lclean_bgacc = 0;
 284     ec->Ltx = ec->Lrx = ec->Lclean = ec->Lclean_bg = 0;
 285     ec->tx_1 = ec->tx_2 = ec->rx_1 = ec->rx_2 = 0;
 286     ec->Lbgn = ec->Lbgn_acc = 0;
 287     ec->Lbgn_upper = 200;
 288     ec->Lbgn_upper_acc = ec->Lbgn_upper << 13;
 289
 290     return  ec;
 291 }
 292 EXPORT_SYMBOL_GPL(oslec_create);
 293 /*- End of function --------------------------------------------------------*/
 294
 295 void oslec_free(struct oslec_state *ec)
 296 {
 297         int i;
 298
 299         fir16_free(&ec->fir_state);
 300         fir16_free(&ec->fir_state_bg);
 301         for (i = 0;  i < 2;  i++)
 302                 kfree(ec->fir_taps16[i]);
 303         kfree(ec->snapshot);
 304         kfree(ec);
 305 }
 306 EXPORT_SYMBOL_GPL(oslec_free);
 307 /*- End of function --------------------------------------------------------*/
 308
 309 void oslec_adaption_mode(struct oslec_state *ec, int adaption_mode)
 310 {
 311     ec->adaption_mode = adaption_mode;
 312 }
 313 EXPORT_SYMBOL_GPL(oslec_adaption_mode);
 314 /*- End of function --------------------------------------------------------*/
 315
 316 void oslec_flush(struct oslec_state *ec)
 317 {
 318     int i;
 319
 320     ec->Ltxacc = ec->Lrxacc = ec->Lcleanacc = ec->Lclean_bgacc = 0;
 321     ec->Ltx = ec->Lrx = ec->Lclean = ec->Lclean_bg = 0;
 322     ec->tx_1 = ec->tx_2 = ec->rx_1 = ec->rx_2 = 0;
 323
 324     ec->Lbgn = ec->Lbgn_acc = 0;
 325     ec->Lbgn_upper = 200;
 326     ec->Lbgn_upper_acc = ec->Lbgn_upper << 13;
 327
 328     ec->nonupdate_dwell = 0;
 329
 330     fir16_flush(&ec->fir_state);
 331     fir16_flush(&ec->fir_state_bg);
 332     ec->fir_state.curr_pos = ec->taps - 1;
 333     ec->fir_state_bg.curr_pos = ec->taps - 1;
 334     for (i = 0;  i < 2;  i++)
 335         memset(ec->fir_taps16[i], 0, ec->taps*sizeof(int16_t));
 336
 337     ec->curr_pos = ec->taps - 1;
 338     ec->Pstates = 0;
 339 }
 340 EXPORT_SYMBOL_GPL(oslec_flush);
 341 /*- End of function --------------------------------------------------------*/
 342
 343 void oslec_snapshot(struct oslec_state *ec) {
 344     memcpy(ec->snapshot, ec->fir_taps16[0], ec->taps*sizeof(int16_t));
 345 }
 346 EXPORT_SYMBOL_GPL(oslec_snapshot);
 347 /*- End of function --------------------------------------------------------*/
 348
 349 /* Dual Path Echo Canceller ------------------------------------------------*/
 350
 351 int16_t oslec_update(struct oslec_state *ec, int16_t tx, int16_t rx)
 352 {
 353     int32_t echo_value;
 354     int clean_bg;
 355     int tmp, tmp1;
 356
 357     /* Input scaling was found be required to prevent problems when tx
 358        starts clipping.  Another possible way to handle this would be the
 359        filter coefficent scaling. */
 360
 361     ec->tx = tx; ec->rx = rx;
 362     tx >>=1;
 363     rx >>=1;
 364
 365     /*
 366        Filter DC, 3dB point is 160Hz (I think), note 32 bit precision required
 367        otherwise values do not track down to 0. Zero at DC, Pole at (1-Beta)
 368        only real axis.  Some chip sets (like Si labs) don't need
 369        this, but something like a $10 X100P card does.  Any DC really slows
 370        down convergence.
 371
 372        Note: removes some low frequency from the signal, this reduces
 373        the speech quality when listening to samples through headphones
 374        but may not be obvious through a telephone handset.
 375
 376        Note that the 3dB frequency in radians is approx Beta, e.g. for
 377        Beta = 2^(-3) = 0.125, 3dB freq is 0.125 rads = 159Hz.
 378     */
 379
 380     if (ec->adaption_mode & ECHO_CAN_USE_RX_HPF) {
 381       tmp = rx << 15;
 382 #if 1
 383         /* Make sure the gain of the HPF is 1.0. This can still saturate a little under
 384            impulse conditions, and it might roll to 32768 and need clipping on sustained peak
 385            level signals. However, the scale of such clipping is small, and the error due to
 386            any saturation should not markedly affect the downstream processing. */
 387         tmp -= (tmp >> 4);
 388 #endif
 389       ec->rx_1 += -(ec->rx_1>>DC_LOG2BETA) + tmp - ec->rx_2;
 390
 391       /* hard limit filter to prevent clipping.  Note that at this stage
 392          rx should be limited to +/- 16383 due to right shift above */
 393       tmp1 = ec->rx_1 >> 15;
 394       if (tmp1 > 16383) tmp1 = 16383;
 395       if (tmp1 < -16383) tmp1 = -16383;
 396       rx = tmp1;
 397       ec->rx_2 = tmp;
 398     }
 399
 400     /* Block average of power in the filter states.  Used for
 401        adaption power calculation. */
 402
 403     {
 404         int new, old;
 405
 406         /* efficient "out with the old and in with the new" algorithm so
 407            we don't have to recalculate over the whole block of
 408            samples. */
 409         new = (int)tx * (int)tx;
 410         old = (int)ec->fir_state.history[ec->fir_state.curr_pos] *
 411               (int)ec->fir_state.history[ec->fir_state.curr_pos];
 412         ec->Pstates += ((new - old) + (1<<ec->log2taps)) >> ec->log2taps;
 413         if (ec->Pstates < 0) ec->Pstates = 0;
 414     }
 415
 416     /* Calculate short term average levels using simple single pole IIRs */
 417
 418     ec->Ltxacc += abs(tx) - ec->Ltx;
 419     ec->Ltx = (ec->Ltxacc + (1<<4)) >> 5;
 420     ec->Lrxacc += abs(rx) - ec->Lrx;
 421     ec->Lrx = (ec->Lrxacc + (1<<4)) >> 5;
 422
 423     /* Foreground filter ---------------------------------------------------*/
 424
 425     ec->fir_state.coeffs = ec->fir_taps16[0];
 426     echo_value = fir16(&ec->fir_state, tx);
 427     ec->clean = rx - echo_value;
 428     ec->Lcleanacc += abs(ec->clean) - ec->Lclean;
 429     ec->Lclean = (ec->Lcleanacc + (1<<4)) >> 5;
 430
 431     /* Background filter ---------------------------------------------------*/
 432
 433     echo_value = fir16(&ec->fir_state_bg, tx);
 434     clean_bg = rx - echo_value;
 435     ec->Lclean_bgacc += abs(clean_bg) - ec->Lclean_bg;
 436     ec->Lclean_bg = (ec->Lclean_bgacc + (1<<4)) >> 5;
 437
 438     /* Background Filter adaption -----------------------------------------*/
 439
 440     /* Almost always adap bg filter, just simple DT and energy
 441        detection to minimise adaption in cases of strong double talk.
 442        However this is not critical for the dual path algorithm.
 443     */
 444     ec->factor = 0;
 445     ec->shift = 0;
 446     if ((ec->nonupdate_dwell == 0)) {
 447         int   P, logP, shift;
 448
 449         /* Determine:
 450
 451            f = Beta * clean_bg_rx/P ------ (1)
 452
 453            where P is the total power in the filter states.
 454
 455            The Boffins have shown that if we obey (1) we converge
 456            quickly and avoid instability.
 457
 458            The correct factor f must be in Q30, as this is the fixed
 459            point format required by the lms_adapt_bg() function,
 460            therefore the scaled version of (1) is:
 461
 462            (2^30) * f  = (2^30) * Beta * clean_bg_rx/P
 463                factor  = (2^30) * Beta * clean_bg_rx/P         ----- (2)
 464
 465            We have chosen Beta = 0.25 by experiment, so:
 466
 467                factor  = (2^30) * (2^-2) * clean_bg_rx/P
 468
 469                                        (30 - 2 - log2(P))
 470                factor  = clean_bg_rx 2                         ----- (3)
 471
 472            To avoid a divide we approximate log2(P) as top_bit(P),
 473            which returns the position of the highest non-zero bit in
 474            P.  This approximation introduces an error as large as a
 475            factor of 2, but the algorithm seems to handle it OK.
 476
 477            Come to think of it a divide may not be a big deal on a
 478            modern DSP, so its probably worth checking out the cycles
 479            for a divide versus a top_bit() implementation.
 480         */
 481
 482         P = MIN_TX_POWER_FOR_ADAPTION + ec->Pstates;
 483         logP = top_bit(P) + ec->log2taps;
 484         shift = 30 - 2 - logP;
 485         ec->shift = shift;
 486
 487         lms_adapt_bg(ec, clean_bg, shift);
 488     }
 489
 490     /* very simple DTD to make sure we dont try and adapt with strong
 491        near end speech */
 492
 493     ec->adapt = 0;
 494     if ((ec->Lrx > MIN_RX_POWER_FOR_ADAPTION) && (ec->Lrx > ec->Ltx))
 495         ec->nonupdate_dwell = DTD_HANGOVER;
 496     if (ec->nonupdate_dwell)
 497         ec->nonupdate_dwell--;
 498
 499     /* Transfer logic ------------------------------------------------------*/
 500
 501     /* These conditions are from the dual path paper [1], I messed with
 502        them a bit to improve performance. */
 503
 504     if ((ec->adaption_mode & ECHO_CAN_USE_ADAPTION) &&
 505         (ec->nonupdate_dwell == 0) &&
 506         (8*ec->Lclean_bg < 7*ec->Lclean) /* (ec->Lclean_bg < 0.875*ec->Lclean) */ &&
 507         (8*ec->Lclean_bg < ec->Ltx)      /* (ec->Lclean_bg < 0.125*ec->Ltx)    */ )
 508     {
 509         if (ec->cond_met == 6) {
 510             /* BG filter has had better results for 6 consecutive samples */
 511             ec->adapt = 1;
 512             memcpy(ec->fir_taps16[0], ec->fir_taps16[1], ec->taps*sizeof(int16_t));
 513         }
 514         else
 515             ec->cond_met++;
 516     }
 517     else
 518         ec->cond_met = 0;
 519
 520     /* Non-Linear Processing ---------------------------------------------------*/
 521
 522     ec->clean_nlp = ec->clean;
 523     if (ec->adaption_mode & ECHO_CAN_USE_NLP)
 524     {
 525         /* Non-linear processor - a fancy way to say "zap small signals, to avoid
 526            residual echo due to (uLaw/ALaw) non-linearity in the channel.". */
 527
 528       if ((16*ec->Lclean < ec->Ltx))
 529       {
 530         /* Our e/c has improved echo by at least 24 dB (each factor of 2 is 6dB,
 531            so 2*2*2*2=16 is the same as 6+6+6+6=24dB) */
 532         if (ec->adaption_mode & ECHO_CAN_USE_CNG)
 533         {
 534             ec->cng_level = ec->Lbgn;
 535
 536             /* Very elementary comfort noise generation.  Just random
 537                numbers rolled off very vaguely Hoth-like.  DR: This
 538                noise doesn't sound quite right to me - I suspect there
 539                are some overlfow issues in the filtering as it's too
 540                "crackly".  TODO: debug this, maybe just play noise at
 541                high level or look at spectrum.
 542             */
 543
 544             ec->cng_rndnum = 1664525U*ec->cng_rndnum + 1013904223U;
 545             ec->cng_filter = ((ec->cng_rndnum & 0xFFFF) - 32768 + 5*ec->cng_filter) >> 3;
 546             ec->clean_nlp = (ec->cng_filter*ec->cng_level*8) >> 14;
 547
 548         }
 549         else if (ec->adaption_mode & ECHO_CAN_USE_CLIP)
 550         {
 551             /* This sounds much better than CNG */
 552             if (ec->clean_nlp > ec->Lbgn)
 553               ec->clean_nlp = ec->Lbgn;
 554             if (ec->clean_nlp < -ec->Lbgn)
 555               ec->clean_nlp = -ec->Lbgn;
 556         }
 557         else
 558         {
 559           /* just mute the residual, doesn't sound very good, used mainly
 560              in G168 tests */
 561           ec->clean_nlp = 0;
 562         }
 563       }
 564       else {
 565           /* Background noise estimator.  I tried a few algorithms
 566              here without much luck.  This very simple one seems to
 567              work best, we just average the level using a slow (1 sec
 568              time const) filter if the current level is less than a
 569              (experimentally derived) constant.  This means we dont
 570              include high level signals like near end speech.  When
 571              combined with CNG or especially CLIP seems to work OK.
 572           */
 573           if (ec->Lclean < 40) {
 574               ec->Lbgn_acc += abs(ec->clean) - ec->Lbgn;
 575               ec->Lbgn = (ec->Lbgn_acc + (1<<11)) >> 12;
 576           }
 577        }
 578     }
 579
 580     /* Roll around the taps buffer */
 581     if (ec->curr_pos <= 0)
 582         ec->curr_pos = ec->taps;
 583     ec->curr_pos--;
 584
 585     if (ec->adaption_mode & ECHO_CAN_DISABLE)
 586       ec->clean_nlp = rx;
 587
 588     /* Output scaled back up again to match input scaling */
 589
 590     return (int16_t) ec->clean_nlp << 1;
 591 }
 592 EXPORT_SYMBOL_GPL(oslec_update);
 593 /*- End of function --------------------------------------------------------*/
 594
 595 /* This function is seperated from the echo canceller is it is usually called
 596    as part of the tx process.  See rx HP (DC blocking) filter above, it's
 597    the same design.
 598
 599    Some soft phones send speech signals with a lot of low frequency
 600    energy, e.g. down to 20Hz.  This can make the hybrid non-linear
 601    which causes the echo canceller to fall over.  This filter can help
 602    by removing any low frequency before it gets to the tx port of the
 603    hybrid.
 604
 605    It can also help by removing and DC in the tx signal.  DC is bad
 606    for LMS algorithms.
 607
 608    This is one of the classic DC removal filters, adjusted to provide sufficient
 609    bass rolloff to meet the above requirement to protect hybrids from things that
 610    upset them. The difference between successive samples produces a lousy HPF, and
 611    then a suitably placed pole flattens things out. The final result is a nicely
 612    rolled off bass end. The filtering is implemented with extended fractional
 613    precision, which noise shapes things, giving very clean DC removal.
 614 */
 615
 616 int16_t oslec_hpf_tx(struct oslec_state *ec, int16_t tx) {
 617     int tmp, tmp1;
 618
 619     if (ec->adaption_mode & ECHO_CAN_USE_TX_HPF) {
 620         tmp = tx << 15;
 621 #if 1
 622         /* Make sure the gain of the HPF is 1.0. The first can still saturate a little under
 623            impulse conditions, and it might roll to 32768 and need clipping on sustained peak
 624            level signals. However, the scale of such clipping is small, and the error due to
 625            any saturation should not markedly affect the downstream processing. */
 626         tmp -= (tmp >> 4);
 627 #endif
 628         ec->tx_1 += -(ec->tx_1>>DC_LOG2BETA) + tmp - ec->tx_2;
 629         tmp1 = ec->tx_1 >> 15;
 630         if (tmp1 > 32767) tmp1 = 32767;
 631         if (tmp1 < -32767) tmp1 = -32767;
 632         tx = tmp1;
 633         ec->tx_2 = tmp;
 634     }
 635
 636     return tx;
 637 }
 638 EXPORT_SYMBOL_GPL(oslec_hpf_tx);
 639
 640 MODULE_LICENSE("GPL");
 641 MODULE_AUTHOR("David Rowe");
 642 MODULE_DESCRIPTION("Open Source Line Echo Canceller");
 643 MODULE_VERSION("0.3.0");