Pastebin launched a little side project called VERYVIRAL.com, check it out ;-) Want more features on Pastebin? Sign Up, it's FREE!
Guest

signal_processing

By: a guest on Feb 11th, 2014  |  syntax: C++  |  size: 5.05 KB  |  views: 37  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. #include <ipp.h>
  2. #include <cmath>
  3. #include <iostream>
  4. #include <fstream>
  5. #include <stdexcept>
  6. #include <memory>
  7.  
  8. // Signal parameter.
  9. Ipp32s signal_frequency [5] = { 1, 5, 10, 200, 500 };   // in Hz.
  10. Ipp32s signal_amplitude = 20;                                                   // in levels.
  11. Ipp32s noise_amplitude = 10;                                                    // in levels.
  12. const Ipp32s sampling_frequency = 44100;                                // in Hz.
  13.  
  14. // Free functor.
  15. struct freeIpps {
  16.         inline void operator()(void *_p)
  17.         {
  18.                 if(!_p) return;
  19.  
  20.                 ippsFree(_p);
  21.         }
  22. };
  23.  
  24. // Type references.
  25. typedef std::unique_ptr<Ipp32f[], freeIpps> Ipp32f_ptr;
  26. typedef std::unique_ptr<Ipp8u[], freeIpps> Ipp8u_ptr;
  27. typedef std::unique_ptr<IppsFFTSpec_R_32f, freeIpps> IppsFFTSpec_R_32f_ptr;
  28.  
  29. // Signal data.
  30. const Ipp32f time_step = 1;
  31. const Ipp32f log_2 = std::log(2);
  32. Ipp32f time_line[sampling_frequency];
  33.  
  34. // Check function.
  35. inline bool ispwr_of_2(Ipp32s x)
  36. {
  37.         return (x && !(x & (x - 1)));
  38. }
  39.  
  40. inline Ipp32s topwr_of_2(Ipp32s x)
  41. {
  42.         if(x > 0)
  43.         {
  44.                 if(x & (x - 1))
  45.                 {
  46.                         --x;
  47.  
  48.                         x |= x >> 1;
  49.                         x |= x >> 2;
  50.                         x |= x >> 4;
  51.                         x |= x >> 8;
  52.                         x |= x >> 16;
  53.  
  54.                         ++x;
  55.                 }
  56.         }      
  57.  
  58.         return  x;
  59. }
  60.  
  61. void init_time()
  62. {
  63.         time_line[0] = time_step;
  64.  
  65.         for(Ipp32u i = 1; i < sampling_frequency; ++i)
  66.                 time_line[i] = time_line[i-1] + time_step;
  67. }
  68.  
  69. Ipp32f_ptr signal_gen(Ipp32s fi, Ipp32s &data_sz)
  70. {
  71.         // The length of the vector transformed by the FFT must be a power of 2.
  72.         data_sz = topwr_of_2(sampling_frequency);
  73.         Ipp32f_ptr data(ippsMalloc_32f(data_sz));
  74.  
  75.         Ipp32f y = IPP_2PI * signal_frequency[fi];
  76.         y /= sampling_frequency;
  77.  
  78.         // Generate sine wave with sampling_frequency size.
  79.         for(Ipp32u i = 0; i < sampling_frequency; ++i) {
  80.                 data[i] = y * time_line[i];
  81.                 data[i] = std::sin(data[i]);
  82.                 data[i] *= signal_amplitude;
  83.         }
  84.  
  85.         Ipp32s zero_begin = sampling_frequency;
  86.         Ipp32s zero_end = data_sz - zero_begin;
  87.  
  88.         // Zero-pad the input signals (add zeros to the end so that at least half of the wave is "blank").
  89.         IppStatus status = ippsZero_32f(data.get() + zero_begin, zero_end);
  90.         if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  91.  
  92.         std::cout << "The sine wave is generated. The length is equal to: " << sampling_frequency
  93.                           << "\nZeroed part of the signal is in the range: [ " << zero_begin << " ] [ " << zero_begin + zero_end << " ]\n";
  94.  
  95.         return data;
  96. }
  97.  
  98. void do_processing()
  99. {
  100.         while(true)
  101.         for(Ipp32u i = 0; i < 5;  i = (i+1 == 5)?(0):(i+1))
  102.         {
  103.                 Ipp32s fftSz = 0; // Length of data with zeros.
  104.                 Ipp32f_ptr pSrc = signal_gen(i, fftSz);
  105.                 Ipp32f_ptr pDst(ippsMalloc_32f(fftSz));
  106.  
  107.                 Ipp32s pSpecSize = 0;
  108.                 Ipp32s pSpecBuffSize = 0;
  109.                 Ipp32s pBuffSize = 0;
  110.                 Ipp32s order = std::log(fftSz)/log_2;
  111.  
  112.                 IppStatus status = ippsFFTGetSize_R_32f(order, IPP_FFT_NODIV_BY_ANY, ippAlgHintFast, &pSpecSize, &pSpecBuffSize, &pBuffSize);
  113.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  114.  
  115.                 Ipp8u_ptr pBuffer(ippsMalloc_8u(pBuffSize));
  116.                 Ipp8u_ptr pSpec(ippsMalloc_8u(pSpecSize));
  117.                 Ipp8u_ptr pSpecBuffer(ippsMalloc_8u(pSpecBuffSize));
  118.  
  119.                 IppsFFTSpec_R_32f *pFFTSpec = nullptr;
  120.  
  121.                 status = ippsFFTInit_R_32f(&pFFTSpec, order, IPP_FFT_NODIV_BY_ANY, ippAlgHintFast, pSpec.get(), pSpecBuffer.get());
  122.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  123.                
  124.                 status = ippsFFTFwd_RToPack_32f(pSrc.get(), pDst.get(), pFFTSpec, pBuffer.get());
  125.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  126.  
  127.                 Ipp32fc *pDstC = (Ipp32fc*)(pDst.get());
  128.                 Ipp32s realSz = fftSz/2;                  // Exclude zero part.        
  129.                 Ipp32s magSz = realSz/2;                  // Convert length to complex length.
  130.  
  131.                 // Extract magniutde from input signal.
  132.                 status = ippsMagnitude_32fc(pDstC, pDst.get(), realSz);
  133.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  134.  
  135.                 Ipp32f lvl = 2.0/Ipp32f(realSz);
  136.                 Ipp32f B = 19.93156857;
  137.  
  138.                 // Get amplitude in levels.
  139.                 status = ippsMulC_32f(pDst.get(), lvl, pDst.get(), magSz);
  140.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  141.  
  142.                 // Convert level to dB.
  143.                 status = ippsLn_32f(pDst.get(), pDst.get(), magSz);
  144.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  145.  
  146.                 status = ippsMulC_32f(pDst.get(), B, pDst.get(), magSz);
  147.                 if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  148.         }      
  149. }
  150.  
  151. int main()
  152. try{
  153.         IppStatus status = ippInit();
  154.         if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  155.  
  156.         std::cout << "The IPP library was successfully initialized.\n";
  157.  
  158.         status = ippSetAffinity(ippAffinityCompactFineHT, 0);
  159.         if(status != ippStsNoErr) throw std::runtime_error(ippGetStatusString(status));
  160.  
  161.         std::cout << "OpenMP succesfully bind all threads to OS processors. \n";
  162.  
  163.         init_time();
  164.         std::cout << "Timeline was initialized.\n";
  165.  
  166.         do_processing();
  167.         std::cout << "Signal processing has completed successfully.\n";
  168.  
  169.         return 0;
  170. }
  171. catch(std::runtime_error &e)
  172. {
  173.         std::cerr << e.what() << '\n';
  174.         return 1;
  175. }
  176. catch(...)
  177. {
  178.         std::cerr << "Unhandled exception.\n";
  179.         return 2;
  180. }