Changeset 4844

Show
Ignore:
Timestamp:
08/27/07 00:05:50 (6 years ago)
Author:
jgaeddert
Message:

fleshing out implementation of P/N sequence generator

Location:
SigProc/trunk/SigProc
Files:
2 modified
2 copied

Legend:

Unmodified
Added
Removed

  • SigProc/trunk/SigProc/Makefile.am

    r4766 r4844  
    44    pll.cpp                         \ 
    55    scaling.cpp                     \ 
     6    sequencing.cpp                  \ 
    67    sources.cpp                     \ 
    78    utility.cpp                     \ 
    … …  
    2930    pll.o                           \ 
    3031    scaling.o                       \ 
     32    sequencing.o                    \ 
    3133    sources.o                       \ 
    3234    utility.o                       \ 
    … …  
    4547        autotest_sources/dot_product_testsuite.h \ 
    4648        autotest_sources/trig_testsuite.h \ 
     49        autotest_sources/sequencing_testsuite.h \ 
    4750        autotest_sources/Fec_conv_testsuite.h 
    4851 

  • SigProc/trunk/SigProc/SigProc.h

    r4766 r4844  
    401401{ 
    402402  public: 
    403     /// default constructor 
    404     PNSequence(); 
     403    /// initializing constructor 
     404    PNSequence(unsigned int _g, unsigned int _a); 
    405405 
    406406    /// destructor 
    407407    ~PNSequence(); 
    408408 
    409   private: 
    410     // g: generator polynomial 
    411     // a: initial polynomial state 
     409    unsigned long m;    ///< shift register length 
     410    unsigned long n;    ///< output sequence length, \f$ n=2^m-1 \f$ 
     411    char * g;           ///< generator polynomial 
     412    char * a;           ///< initial polynomial state 
     413    char * s;           ///< output sequence 
    412414}; 
    413415 

  • SigProc/trunk/SigProc/autotest_sources/sequencing_testsuite.h

    r4766 r4844  
    1 #ifndef __TRIG_TEST_H__ 
    2 #define __TRIG_TEST_H__ 
     1#ifndef __SEQUENCING_TEST_H__ 
     2#define __SEQUENCING_TEST_H__ 
    33 
    44#include <cxxtest/TestSuite.h> 
    … …  
    99// 
    1010 
    11 class arctan_Testsuite : public CxxTest::TestSuite 
     11class PNSequence_Testsuite : public CxxTest::TestSuite 
    1212{ 
    1313public: 
    1414 
    15     // 
    16     // arctan( float&, float&, float&) 
    17     // 
    18      
    19     void test_float_float_float_limits() { 
    20         float x, y, theta; 
     15    void test_constructor_01() { 
     16        unsigned long g(0x0007); 
     17        unsigned long a(0x0001); 
     18        SigProc::PNSequence pn(g, a); 
     19        TS_ASSERT_EQUALS(pn.m, 2); 
     20        TS_ASSERT_EQUALS(pn.n, 3); 
    2121 
    22         x = 0.0f; 
    23         y = 1.0f; 
    24         arctan(x, y, theta); 
    25         TS_ASSERT_DELTA(theta, PI/2, 1e-9); 
    26  
    27         x = 0.0f; 
    28         y = -1.0f; 
    29         arctan(x, y, theta); 
    30         TS_ASSERT_DELTA(theta, 3*PI/2, 1e-9); 
    31  
    32         x = 0.0f; 
    33         y = 0.0f; 
    34         arctan(x, y, theta); 
    35         TS_ASSERT_DELTA(theta, PI/2, 1e-9); 
    36  
    37         x = 1.0; 
    38         y = 0.0f; 
    39         arctan(x, y, theta); 
    40         TS_ASSERT_DELTA(theta, 0, 1e-9); 
    41  
     22        // TODO: Generate output sequence 
    4223    } 
    4324 
     25    void test_constructor_02() { 
     26        unsigned long g(0x0043); 
     27        unsigned long a(0x0001); 
     28        SigProc::PNSequence pn(g, a); 
     29        TS_ASSERT_EQUALS(pn.m, 6); 
     30        TS_ASSERT_EQUALS(pn.n, 63); 
    4431 
    45     void test_float_float_float_quadrants() { 
    46         float x, y, theta; 
    47  
    48         // Q1 
    49         x = 1.0f; 
    50         y = 1.0f; 
    51         arctan(x, y, theta); 
    52         TS_ASSERT_DELTA(theta, PI/4, 1e-9); 
    53  
    54         // Q2 
    55         x = -1.0f; 
    56         y = 1.0f; 
    57         arctan(x, y, theta); 
    58         TS_ASSERT_DELTA(theta, 3*PI/4, 1e-9); 
    59  
    60         // Q3 
    61         x = -1.0f; 
    62         y = -1.0f; 
    63         arctan(x, y, theta); 
    64         TS_ASSERT_DELTA(theta, 5*PI/4, 1e-9); 
    65  
    66         // Q4 
    67         x = 1.0f; 
    68         y = -1.0f; 
    69         arctan(x, y, theta); 
    70         TS_ASSERT_DELTA(theta, 7*PI/4, 1e-9); 
    71  
     32        // TODO: Generate output sequence 
    7233    } 
    73  
    7434 
    7535}; 
    7636 
    7737 
    78 #endif // __TRIG_TEST_H__ 
     38#endif // __SEQUENCING_TEST_H__ 
    7939 

  • SigProc/trunk/SigProc/sequencing.cpp

    r4766 r4844  
    2121****************************************************************************/ 
    2222 
     23/*! \file This file describes the bit sequencing classes and functions 
     24 * 
     25 */ 
     26 
    2327#include "SigProc.h" 
    2428 
    … …  
    2731//----------------------------------------------------------------------------- 
    2832// 
    29 // Design root raised-cosine filter 
    30 // 
    31 //----------------------------------------------------------------------------- 
    32 /** \brief Calculate square-root raised-cosine filter coefficients 
    33  
    34   DesignRRCFilter calculates the coefficients for a square-root raised-cosine 
    35   (RRC) finite impulse response (FIR) filter commonly used in digital 
    36   communications.  The input parameters are as follows 
    37  
    38     - \f$k\f$ : samples per symbol 
    39     - \f$m\f$ : sample delay 
    40     - \f$\beta\f$ : excess bandwidth (rolloff) factor 
    41  
    42   The function returns a pointer to the filter coefficients as well as an 
    43   integer value describing the length of the filter.  The length of the 
    44   filter is always 
    45  
    46   \f[ 
    47     h_{len} = 2 k m + 1 
    48   \f] 
    49  
    50   The filter coefficients themselves are derived from the following equation 
    51  
    52   \f[  
    53     h\left[z\right] = 
    54       4\beta \frac{ \cos\left[(1+\beta)\pi z\right] + 
    55                     \sin\left[(1+\beta)\pi z\right] / (4\beta z) } 
    56                   { \pi \sqrt{T}\left[ 1-16\beta^2z^2\right] } 
    57   \f] 
    58  
    59   where \f$z=n/k-m\f$, and \f$T=1\f$ for most cases. 
    60  
    61   The function compensates for the two cases where \f$h[n]\f$ might be 
    62   undefined in the above equation, viz. 
    63  
    64   \f[ 
    65     \mathop {\lim }\limits_{z \to 0 } h(z) = 
    66       \frac{ 4\beta \left[ 1 + \frac{1-\beta\pi }{ 4\beta } \right] } 
    67            { \pi\sqrt{T}\left( 1-16\beta^2 z^2 \right) } 
    68   \f] 
    69  
    70   and 
    71  
    72   \f[ 
    73     \mathop {\lim }\limits_{z \to \pm \frac{1}{4\beta} } h(z) = 
    74         \frac{(1+\beta)}{2\pi}\sin\left[\frac{(1+\beta)\pi}{4\beta}\right] 
    75       - \frac{(1-\beta)}{2}\cos\left[\frac{(1-\beta)\pi}{4\beta}\right] 
    76       + \frac{2\beta}{\pi}\sin\left[\frac{(1-\beta)\pi}{4\beta}\right] 
    77  
    78   \f] 
    79  
    80   \param[in]  k         samples per symbol 
    81   \param[in]  m         symbol delay 
    82   \param[in]  beta      excess bandwidth/rolloff factor ( 0 < beta < 1 ) 
    83   \param[out] h         pointer to filter coefficients 
    84   \param[out] h_len     length of filter, h_len = 2*m*k+1 
    85  
    86  */ 
    87  
    88 void DesignRRCFilter( 
    89   unsigned int k,      // samples per symbol 
    90   unsigned int m,      // delay 
    91   float beta,          // rolloff factor ( 0 < beta <= 1 ) 
    92   float *& h,          // pointer to filter coefficients 
    93   unsigned int & h_len // length of filter (len = 2*m*k+1) 
    94 ) 
    95 { 
    96     h_len = 0; 
    97  
    98     if ( k < 1 ) { 
    99         std::cerr << "ERROR: SigProc::DesignRRCFilter: k must be greater than 0" 
    100                   << std::endl; 
    101         throw 0; 
    102     } else if ( m < 1 ) { 
    103         std::cerr << "ERROR: SigProc::DesignRRCFilter: m must be greater than 0" 
    104                   << std::endl; 
    105         throw 0; 
    106     } else if ( (beta < 0.0f) || (beta > 1.0f) ) { 
    107         std::cerr << "ERROR: SigProc::DesignRRCFilter: beta must be in [0,1]" 
    108                   << std::endl; 
    109         throw 0; 
    110     } else; 
    111  
    112     unsigned int n; 
    113     float z, t1, t2, t3, t4, T(1.0f); 
    114     float pi(3.14159265358979f); 
    115  
    116     h_len = 2*k*m+1; 
    117     h = new float[h_len]; 
    118  
    119     // Calculate filter coefficients 
    120     for (n=0; n<h_len; n++) { 
    121  
    122         z = float(n)/float(k)-float(m); 
    123         t1 = cosf((1+beta)*pi*z); 
    124         t2 = sinf((1-beta)*pi*z); 
    125  
    126         // Check for special condition where z equals zero 
    127         if ( n == k*m ) { 
    128             t4 = 4*beta/(pi*sqrtf(T)*(1-(16*beta*beta*z*z))); 
    129             h[n] = t4*( 1 + (1-beta)*pi/(4*beta) ); 
    130         } else { 
    131             t3 = 1/((4*beta*z)); 
    132  
    133             float g = 1-16*beta*beta*z*z; 
    134             g *= g; 
    135  
    136             // Check for special condition where 16*beta^2*z^2 equals 1 
    137             if ( g < 1e-3 ) { 
    138                 float g1, g2, g3, g4; 
    139                 g1 = -(1+beta)*pi*sin((1+beta)*pi/(4*beta)); 
    140                 g2 = cos((1-beta)*pi/(4*beta))*(1-beta)*pi; 
    141                 g3 = -sin((1-beta)*pi/(4*beta))*4*beta; 
    142                 g4 = -2*pi; 
    143  
    144                 h[n] = (g1+g2+g3)/g4; 
    145             } else { 
    146                 t4 = 4*beta/(pi*sqrtf(T)*(1-(16*beta*beta*z*z))); 
    147                 h[n] = t4*( t1 + (t2*t3) ); 
    148             } 
    149         } 
    150     } 
    151 } 
    152  
    153 //----------------------------------------------------------------------------- 
    154 // 
    155 // Design Gaussian filter 
    156 // 
    157 //----------------------------------------------------------------------------- 
    158 /** \brief Calculate Gaussian filter coefficients 
    159  
    160   DesignRRCFilter calculates the coefficients for a square-root raised-cosine 
    161   (RRC) finite impulse response (FIR) filter commonly used in digital 
    162   communications.  The input parameters are as follows 
    163  
    164     - \f$k\f$ : samples per symbol 
    165     - \f$m\f$ : sample delay 
    166     - \f$\beta\f$ : excess bandwidth factor 
    167  
    168   The function returns a pointer to the filter coefficients as well as an 
    169   integer value describing the length of the filter.  The length of the 
    170   filter is always 
    171  
    172   \f[ 
    173     h_{len} = 2 k m + 1 
    174   \f] 
    175  
    176   The filter coefficients themselves are derived from the following equation 
    177  
    178   \f[  
    179     h\left[z\right] = \frac{1}{\sqrt[4]{2\pi\sigma}} 
    180                       e^{\left(-z^2/{4\sigma^2}\right)} 
    181   \f] 
    182  
    183   where \f$z=n/k-m\f$ and \f$\sigma=\beta\f$. 
    184  
    185   \param[in]  k         samples per symbol 
    186   \param[in]  m         symbol delay 
    187   \param[in]  beta      excess bandwidth/rolloff factor ( 0 < beta < 1 ) 
    188   \param[out] h         pointer to filter coefficients 
    189   \param[out] h_len     length of filter, h_len = 2*m*k+1 
    190  
    191  */ 
    192  
    193 void DesignGaussianFilter( 
    194   unsigned int k,      // samples per symbol 
    195   unsigned int m,      // delay 
    196   float beta,          // rolloff factor ( 0 < beta ) 
    197   float *& h,          // pointer to filter coefficients 
    198   unsigned int & h_len // length of filter (len = 2*m*k+1) 
    199 ) { 
    200     h_len = 0; 
    201  
    202     if ( k < 1 ) { 
    203         std::cerr << "ERROR: SigProc::DesignGaussianFilter: k must be greater than 0" 
    204                   << std::endl; 
    205         throw 0; 
    206     } else if ( m < 1 ) { 
    207         std::cerr << "ERROR: SigProc::DesignGaussianFilter: m must be greater than 0" 
    208                   << std::endl; 
    209         throw 0; 
    210     } else if ( beta < 0.0f) { 
    211         std::cerr << "ERROR: SigProc::DesignGaussianFilter: beta must be greater than 0" 
    212                   << std::endl; 
    213         throw 0; 
    214     } else; 
    215  
    216     h_len = 2*k*m + 1; 
    217     h = new float[h_len]; 
    218     float sigma(beta);  ///\todo determine how to calculate sigma from beta 
    219     float pi(3.14159265358979f); 
    220     float g1( 2*sigma*sigma ); 
    221     float g2( 1/sqrtf(g1*pi) ); 
    222     float z; 
    223  
    224     for (unsigned int n=0; n<h_len; n++) { 
    225         // Check for special condition where z equals zero 
    226         if ( n == k*m ) { 
    227             h[n] = g2; 
    228         } else { 
    229             z = float(n)/float(k)-float(m); 
    230             h[n] = expf(-z*z/g1)*g2; 
    231         } 
    232     } 
    233 } 
    234  
    235 //----------------------------------------------------------------------------- 
    236 // 
    237 // FIR polyphase filter bank 
     33// P/N Sequence 
    23834// 
    23935//----------------------------------------------------------------------------- 
    24036 
    241 // Initializing constructor 
    242 FIRPolyphaseFilterBank::FIRPolyphaseFilterBank( 
    243         char * _type,       // type of filter 
    244         unsigned int _k,    // samples per symbol 
    245         unsigned int _m,    // delay 
    246         float _beta,        // excess bandwidth factor 
    247         unsigned int _Npfb  // number of filters 
    248         ) 
    249 { 
    250     k = _k; 
    251     m = _m; 
    252     beta = _beta; 
    253     Npfb = _Npfb; 
    254  
    255     H = NULL; 
    256  
    257     if ( strcmp(_type,"rrcos")==0 ) { 
    258         // Square-root raised-cosine filter 
    259         CalculateRRCFilterCoefficients(); 
    260         TransposeCoefficientMatrix(); 
    261     } else if ( strcmp(_type,"drrcos")==0 ) { 
    262         // Derivative square-root raised-cosine filter 
    263         CalculateRRCFilterCoefficients(); 
    264         CalculateDerivativeFilterCoefficients(); 
    265         TransposeCoefficientMatrix(); 
    266     } else if ( strcmp(_type,"gaussian")==0 ) { 
    267         // Gaussian filter 
    268         CalculateGaussianFilterCoefficients(); 
    269         TransposeCoefficientMatrix(); 
    270     } else if ( strcmp(_type,"dgaussian")==0 ) { 
    271         // Derivative gaussian filter 
    272         CalculateGaussianFilterCoefficients(); 
    273         CalculateDerivativeFilterCoefficients(); 
    274         TransposeCoefficientMatrix(); 
    275     } else { 
    276         std::cerr << "ERROR: FIRPolyphaseFilterBank: unknown filter type : " << _type << std::endl; 
    277         throw 0; 
     37PNSequence::PNSequence(unsigned int _g, unsigned int _a) { 
     38    unsigned int i; 
     39     
     40    // extract shift register length from generator polynomial 
     41    unsigned long g_tmp(_g); 
     42    m = 0; 
     43    // this loop effectively counts the placement of the MSB in _g 
     44    for (i=0; i<sizeof(unsigned long)*8; i++) { 
     45        if ( g_tmp & 0x0001 ) 
     46            m = i; 
     47        g_tmp >>= 1; 
    27848    } 
    279  
    280     // At this point, H should be initialized and h_len should store 
    281     // the correct length of each filter in the bank. 
    282     v.SetBufferSize(h_len); 
    283  
    284     // memset input buffer to zeros 
    285     ResetBuffer(); 
    286 } 
    287  
    288 // Initializing constructor 
    289 FIRPolyphaseFilterBank::FIRPolyphaseFilterBank( 
    290         float * _H,         // filter bank coefficients 
    291         unsigned int _h_len,// length of each filter 
    292         unsigned int _Npfb  // number of filters 
    293         ) 
    294 { 
    295     /// \todo perform check on input? 
    296  
    297     h_len = _h_len; 
    298     Npfb = _Npfb; 
    299  
    300     // perform deep copy of memory 
    301     H = new float[Npfb*h_len]; 
    302     for (unsigned int i=0; i<(Npfb*h_len); i++) 
    303         H[i] = _H[i]; 
    304  
    305     v.SetBufferSize(h_len); 
    306  
    307     // memset input buffer to zeros 
    308     ResetBuffer(); 
    309 } 
    310  
    311 // destructor 
    312 FIRPolyphaseFilterBank::~FIRPolyphaseFilterBank() 
    313 { 
    314     // delete filter bank coefficients matrix 
    315     if ( H != NULL ) 
    316         delete [] H; 
     49     
     50    n = 1; 
     51    for (i=0; i<m; i++) 
     52        n <<= 1; 
     53    n--; 
     54     
     55    // generating polynomial 
     56    g = new char[m]; 
     57    for (i=0; i<m; i++) { 
     58        g[m-i-1] = ( _g & 0x0001 ) ? 1 : 0; 
     59        _g >>= 1; 
     60    } 
     61     
     62    // initial polynomial state 
     63    a = new char[m]; 
     64    for (i=0; i<m; i++) { 
     65        a[m-i-1] = ( _a & 0x0001 ) ? 1 : 0; 
     66        _a >>= 1; 
     67    } 
     68     
     69    // output 
     70    s = new char[n]; 
     71    memset(s, 0x00, n); 
    31772 
    31873} 
    31974 
    320 // push input value into buffer 
    321 void FIRPolyphaseFilterBank::PushInput(short _x) 
    322 { 
    323     // Add _x to input buffer 
    324     v.Push( _x ); 
    325 } 
    326  
    327 // compute filter output from current buffer state using specific filter 
    328 void FIRPolyphaseFilterBank::ComputeOutput( 
    329         short &y,           // output sample 
    330         unsigned int _b     // filter bank index 
    331         ) 
    332 { 
    333     // Ensure the requested filter bank index does not exceed the number 
    334     // of filters actually in the bank 
    335     if ( _b >= Npfb ) 
    336     { 
    337         std::cerr << "ERROR: SigProc::FIRPolyphaseFitlerBank::ComputeOutput, " 
    338                   << " index exceeds filter bank size (" 
    339                   << _b << " > " << Npfb << ")" << std::endl; 
    340         throw 0; 
    341     } 
    342  
    343     // Set B to memory block in filter bank array 
    344     unsigned int B; 
    345     B = _b*h_len; 
    346  
    347     // Compute dot product 
    348     dot_product( H+B, v.GetHeadPtr(), h_len, y); 
    349  
    350 } 
    351  
    352 // compute filter output from current buffer state using specific filter 
    353 void FIRPolyphaseFilterBank::ComputeOutput( 
    354         short &y,           // output sample 
    355         unsigned int _b,    // filter bank index 
    356         short *_v           // external memory buffer array 
    357         ) 
    358 { 
    359     // Ensure the requested filter bank index does not exceed the number 
    360     // of filters actually in the bank 
    361     if ( _b >= Npfb ) 
    362     { 
    363         std::cerr << "ERROR: SigProc::FIRPolyphaseFitlerBank::ComputeOutput, " 
    364                   << " index exceeds filter bank size (" 
    365                   << _b << " > " << Npfb << ")" << std::endl; 
    366         throw 0; 
    367     } 
    368  
    369     // Set B to memory block in filter bank array 
    370     unsigned int B; 
    371     B = _b*h_len; 
    372  
    373     // Compute dot product 
    374     dot_product( H+B, _v, h_len, y); 
    375  
    376 } 
    377  
    378 // Reset filter buffer 
    379 void FIRPolyphaseFilterBank::ResetBuffer() 
    380 { 
    381     v.Release(); 
    382  
    383     // Load buffer with zeros 
    384     for (unsigned int i=0; i<h_len; i++) 
    385         v.Push( 0 ); 
    386 } 
    387  
    388 // Print filter buffer 
    389 void FIRPolyphaseFilterBank::PrintBuffer() 
    390 { 
    391     v.Print(); 
    392 } 
    393  
    394 // Print filter bank coefficients to the screen 
    395 void FIRPolyphaseFilterBank::PrintFilterBankCoefficients() 
    396 { 
    397     if ( H == NULL ) 
    398     { 
    399         std::cout << "ERROR: SigProc::FIRPolyphaseFilterBank: " 
    400                   << "cannot print filter coefficients (matrix empty)" 
    401                   << std::endl; 
    402         return; 
    403     } 
    404  
    405     unsigned int r, c, B; 
    406     std::cout << "H = [" << std::endl; 
    407     for (r=0; r<Npfb; r++ ) { 
    408         B = r*h_len; 
    409         printf("  "); 
    410  
    411         for (c=0; c<h_len; c++) { 
    412             printf("%0.3f ", H[B + c]); 
    413         } 
    414         printf(";\n"); 
    415     } 
    416     std::cout << "   ]" << std::endl; 
    417  
    418 } 
    419  
    420 // 
    421 void FIRPolyphaseFilterBank::TransposeCoefficientMatrix() 
    422 { 
    423     // create temporary pointer 
    424     float * H_tmp; 
    425     H_tmp = H; 
    426  
    427     // allocate new memory for filter bank coefficients 
    428     H = new float[h_len*Npfb]; 
    429  
    430     // reshape matrix 
    431     unsigned int i(0), r, c; 
    432     for (r=0; r<Npfb; r++) { 
    433         for (c=0; c<h_len; c++) 
    434             H[i++] = H_tmp[c*Npfb + r]; 
    435     } 
    436  
    437     delete [] H_tmp; 
    438 } 
    439  
    440  
    441 // Calculate root raised-cosine coefficients 
    442 void FIRPolyphaseFilterBank::CalculateRRCFilterCoefficients() 
    443 { 
    444     if ( H != NULL ) 
    445         delete [] H; 
    446  
    447     // create over-sampled pulse 
    448     DesignRRCFilter(Npfb*k, m, beta, H, h_len); 
    449  
    450     // Apply scaling factor to filter coefficients 
    451     float h_scale = float(k); 
    452     for (unsigned int i=0; i<h_len; i++) 
    453         H[i] /= h_scale; 
    454  
    455     //  
    456     h_len = ( h_len - 1 ) / Npfb; 
    457  
    458 } 
    459  
    460 // 
    461 void FIRPolyphaseFilterBank::CalculateGaussianFilterCoefficients() 
    462 { 
    463     if ( H != NULL ) 
    464         delete [] H; 
    465  
    466 } 
    467  
    468 // Calculate derivative filter coefficients 
    469 void FIRPolyphaseFilterBank::CalculateDerivativeFilterCoefficients() 
    470 { 
    471     unsigned int N = h_len*Npfb; 
    472  
    473     float * dH = new float[N]; 
    474  
    475     for (unsigned int i=0; i<N; i++) { 
    476         if ( i==0 ) { 
    477             dH[0] = H[1] - H[N-1]; 
    478         } else if ( i==N-1 ) { 
    479             dH[N-1] = H[0] - H[N-2]; 
    480         } else { 
    481             dH[i] = H[i+1] - H[i-1]; 
    482         } 
    483         dH[i] /= 2.0f; 
    484     } 
    485  
    486     delete [] H; 
    487     H = dH; 
    488 } 
    489  
    490  
    491 iir_filter::iir_filter(float a_coeff[], unsigned int len_a, float b_coeff[], unsigned int len_b) : len_A(len_a), len_B(len_b), next_v(0) 
    492 { 
    493  
    494     A = new float[len_a]; 
    495     B = new float[len_b]; 
    496  
    497     for (unsigned int i = 0; i < len_a; ++i) 
    498         A[i] = a_coeff[i]; 
    499  
    500     for (unsigned int i = 0; i < len_b; ++i) 
    501         B[i] = b_coeff[i]; 
    502  
    503     if (len_A > len_B) 
    504         len_v = len_A; 
    505     else 
    506         len_v = len_B; 
    507  
    508     v = new float[len_v]; 
    509  
    510     for (unsigned int i = 0; i < len_v; ++i) 
    511         v[i] = 0; 
    512  
    513 } 
    514  
    515 void iir_filter::ResetBuffer() 
    516 { 
    517     for (unsigned int i = 0; i < len_v; ++i) 
    518         v[i] = 0; 
    519 } 
    520  
    521 void iir_filter::do_work(short x, short &y) 
    522 { 
    523     // calculate new v[n] 
    524     int v_idx = next_v; 
    525  
    526     v[next_v] = x; 
    527  
    528     for (unsigned int i = 1; i < len_A; ++i) { 
    529  
    530         --v_idx; 
    531         if (v_idx < 0) 
    532             v_idx += len_v; 
    533  
    534         v[next_v] -= (short) (A[i] * v[v_idx]); 
    535            
    536     } 
    537  
    538     // Now calculate the output value 
    539     v_idx = next_v; 
    540     float out = 0; 
    541  
    542     for (unsigned int i = 0; i < len_B; ++i) { 
    543         out += (B[i] * v[v_idx]); 
    544  
    545         --v_idx; 
    546         if (v_idx < 0) 
    547             v_idx += len_v; 
    548     } 
    549  
    550     next_v = (++next_v) % len_v; 
    551  
    552     if (out < SHRT_MIN) 
    553         y = SHRT_MIN; 
    554     else if (out > SHRT_MAX) 
    555         y = SHRT_MAX; 
    556     else 
    557         y = (short) out; 
    558 } 
    559  
    560  
    561 fir_filter::fir_filter(float a_coeff[], unsigned int len_a) : len_A(len_a), len_v(len_a),next_v(0) 
    562  { 
    563 #ifdef FPM 
    564     A = new mad_fixed_t[len_a]; 
    565     v = new mad_fixed_t[len_v]; 
    566 #else 
    567     A = new float[len_a]; 
    568     v = new short[len_v]; 
    569 #endif 
    570  
    571     for (unsigned int i = 0; i < len_a; ++i) 
    572 #ifdef FPM 
    573         A[i] = mad_f_tofixed(a_coeff[i]); 
    574 #else 
    575         A[i] = a_coeff[i]; 
    576 #endif 
    577  
    578  
    579     for (unsigned int i = 0; i < len_v; ++i) 
    580         v[i] = 0; 
    581  
    582  
    583 } 
    584  
    585 void fir_filter::do_work(bool run_filter, short x, short &y) 
    586 { 
    587     // calculate new v[n] 
    588     int v_idx = next_v; 
    589  
    590     v[next_v] = x; 
    591     next_v = (++next_v) % len_v; 
    592  
    593     if (!run_filter) { // Do not create an output value 
    594         y = 0; 
    595         return; 
    596     } 
    597  
    598     // Now calculate the output value 
    599 #ifdef FPM 
    600     mad_fixed_t out(0); 
    601 #else 
    602     float out(0.0); 
    603 #endif 
    604  
    605     for (unsigned int i = 0; i < len_A; ++i) { 
    606         if (v[v_idx]) { 
    607 #ifdef FPM 
    608             out = mad_f_add(out, mad_f_mul(A[i], v[v_idx])); 
    609 #else 
    610             out += (A[i] * v[v_idx]); 
    611 #endif 
    612         } 
    613  
    614         --v_idx; 
    615         if (v_idx < 0) 
    616             v_idx += len_v; 
    617     } 
    618  
    619 #ifdef FPM 
    620         y = (short) out; 
    621 #else 
    622     if (out < SHRT_MIN) 
    623         y = SHRT_MIN; 
    624     else if (out > SHRT_MAX) 
    625         y = SHRT_MAX; 
    626     else 
    627         y = (short) out; 
    628 #endif 
    629 } 
    630  
    631 void fir_filter::reset() 
    632 { 
    633  
    634     for (unsigned int i = 0; i < len_v; ++i) 
    635         v[i] = 0; 
    636  
    637     next_v = 0; 
     75PNSequence::~PNSequence() { 
     76    delete [] g; 
     77    delete [] a; 
     78    delete [] s; 
    63879} 
    63980 
    64081} // namespace SigProc 
     82