C++ FFT::fwd方法代码示例

void test_scalar_generic(int nfft)
{
    typedef typename FFT<T>::Complex Complex;
    typedef typename FFT<T>::Scalar Scalar;
    typedef typename VectorType<Container, Scalar>::type ScalarVector;
    typedef typename VectorType<Container, Complex>::type ComplexVector;

    FFT<T>        fft;
    ScalarVector  tbuf(nfft);
    ComplexVector freqBuf;

    for (int k = 0; k < nfft; ++k)
        tbuf[k] = (T)(rand() / (double)RAND_MAX - .5);

    // make sure it DOESN'T give the right full spectrum answer
    // if we've asked for half-spectrum
    fft.SetFlag(fft.HalfSpectrum);
    fft.fwd(freqBuf, tbuf);
    VERIFY((size_t)freqBuf.size() == (size_t)((nfft >> 1) + 1));
    VERIFY(fft_rmse(freqBuf, tbuf) < test_precision<T>());  // gross check

    fft.ClearFlag(fft.HalfSpectrum);
    fft.fwd(freqBuf, tbuf);
    VERIFY((size_t)freqBuf.size() == (size_t)nfft);
    VERIFY(fft_rmse(freqBuf, tbuf) < test_precision<T>());  // gross check

    if (nfft & 1)
        return; // odd FFTs get the wrong size inverse FFT

    ScalarVector tbuf2;
    fft.inv(tbuf2, freqBuf);
    VERIFY(dif_rmse(tbuf, tbuf2) < test_precision<T>());  // gross check


    // verify that the Unscaled flag takes effect
    ScalarVector tbuf3;
    fft.SetFlag(fft.Unscaled);

    fft.inv(tbuf3, freqBuf);

    for (int k = 0; k < nfft; ++k)
        tbuf3[k] *= T(1. / nfft);


    // for (size_t i=0;i<(size_t) tbuf.size();++i)
    //    cout << "freqBuf=" << freqBuf[i] << " in2=" << tbuf3[i] << " -  in=" << tbuf[i] << " => " << (tbuf3[i] - tbuf[i] ) <<  endl;

    VERIFY(dif_rmse(tbuf, tbuf3) < test_precision<T>());  // gross check

    // verify that ClearFlag works
    fft.ClearFlag(fft.Unscaled);
    fft.inv(tbuf2, freqBuf);
    VERIFY(dif_rmse(tbuf, tbuf2) < test_precision<T>());  // gross check
}

void bench(int nfft,bool fwd,bool unscaled=false, bool halfspec=false)
{
    typedef typename NumTraits<T>::Real Scalar;
    typedef typename std::complex<Scalar> Complex;
    int nits = NDATA/nfft;
    vector<T> inbuf(nfft);
    vector<Complex > outbuf(nfft);
    FFT< Scalar > fft;

    if (unscaled) {
        fft.SetFlag(fft.Unscaled);
        cout << "unscaled ";
    }
    if (halfspec) {
        fft.SetFlag(fft.HalfSpectrum);
        cout << "halfspec ";
    }


    std::fill(inbuf.begin(),inbuf.end(),0);
    fft.fwd( outbuf , inbuf);

    BenchTimer timer;
    timer.reset();
    for (int k=0;k<8;++k) {
        timer.start();
        if (fwd)
            for(int i = 0; i < nits; i++)
                fft.fwd( outbuf , inbuf);
        else
            for(int i = 0; i < nits; i++)
                fft.inv(inbuf,outbuf);
        timer.stop();
    }

    cout << nameof<Scalar>() << " ";
    double mflops = 5.*nfft*log2((double)nfft) / (1e6 * timer.value() / (double)nits );
    if ( NumTraits<T>::IsComplex ) {
        cout << "complex";
    }else{
        cout << "real   ";
        mflops /= 2;
    }


    if (fwd)
        cout << " fwd";
    else
        cout << " inv";

    cout << " NFFT=" << nfft << "  " << (double(1e-6*nfft*nits)/timer.value()) << " MS/s  " << mflops << "MFLOPS\n";
}

void test_return_by_value(int len)
{
    VectorXf in;
    VectorXf in1;
    in.setRandom( len );
    VectorXcf out1,out2;
    FFT<float> fft;

    fft.SetFlag(fft.HalfSpectrum );

    fft.fwd(out1,in);
    out2 = fft.fwd(in);
    VERIFY( (out1-out2).norm() < test_precision<float>() );
    in1 = fft.inv(out1);
    VERIFY( (in1-in).norm() < test_precision<float>() );
}

void test_complex_generic(int nfft)
{
    typedef typename FFT<T>::Complex Complex;
    typedef typename VectorType<Container,Complex>::type ComplexVector;

    FFT<T> fft;

    ComplexVector inbuf(nfft);
    ComplexVector outbuf;
    ComplexVector buf3;
    for (int k=0;k<nfft;++k)
        inbuf[k]= Complex( (T)(rand()/(double)RAND_MAX - .5), (T)(rand()/(double)RAND_MAX - .5) );
    fft.fwd( outbuf , inbuf);

    VERIFY( fft_rmse(outbuf,inbuf) < test_precision<T>()  );// gross check
    fft.inv( buf3 , outbuf);

    VERIFY( dif_rmse(inbuf,buf3) < test_precision<T>()  );// gross check

    // verify that the Unscaled flag takes effect
    ComplexVector buf4;
    fft.SetFlag(fft.Unscaled);
    fft.inv( buf4 , outbuf);
    for (int k=0;k<nfft;++k)
        buf4[k] *= T(1./nfft);
    VERIFY( dif_rmse(inbuf,buf4) < test_precision<T>()  );// gross check

    // verify that ClearFlag works
    fft.ClearFlag(fft.Unscaled);
    fft.inv( buf3 , outbuf);
    VERIFY( dif_rmse(inbuf,buf3) < test_precision<T>()  );// gross check
}

void fft_Z(MultiArray<complex<float>, 3> &a) {
    FFT<float> fft;
    int nz = a.dims()[2];
    for (int x = 0; x < a.dims()[0]; x++) {
        for (int y = 0; y < a.dims()[1]; y++) {
            VectorXcf fft_in = a.slice<1>({x,y,0},{0,0,nz}).asArray();
            VectorXcf fft_out(nz);
            fft.fwd(fft_out, fft_in);
            a.slice<1>({x,y,0},{0,0,nz}).asArray() = fft_out;
        }
    }
}

void fft_Y(MultiArray<complex<float>, 3> &a) {
    FFT<float> fft;
    int ny = a.dims()[1];
    for (int z = 0; z < a.dims()[2]; z++) {
        for (int x = 0; x < a.dims()[0]; x++) {
            VectorXcf fft_in = a.slice<1>({x,0,z},{0,ny,0}).asArray();
            VectorXcf fft_out(ny);
            fft.fwd(fft_out, fft_in);
            a.slice<1>({x,0,z},{0,ny,0}).asArray() = fft_out;
        }
    }
}

void fft_X(MultiArray<complex<float>, 3> &a) {
    FFT<float> fft;
    int nx = a.dims()[0];
    for (int z = 0; z < a.dims()[2]; z++) {
        for (int y = 0; y < a.dims()[1]; y++) {
            VectorXcf fft_in = a.slice<1>({0,y,z},{nx,0,0}).asArray();
            VectorXcf fft_out(nx);
            fft.fwd(fft_out, fft_in);
            a.slice<1>({0,y,z},{nx,0,0}).asArray() = fft_out;
        }
    }
}

void fwd_inv(size_t nfft)
{
    typedef typename NumTraits<T_freq>::Real Scalar;
    vector<T_time> timebuf(nfft);
    RandomFill(timebuf);

    vector<T_freq> freqbuf;
    static FFT<Scalar> fft;
    fft.fwd(freqbuf,timebuf);

    vector<T_time> timebuf2;
    fft.inv(timebuf2,freqbuf);

    long double rmse = mag2(timebuf - timebuf2) / mag2(timebuf);
    cout << "roundtrip rmse: " << rmse << endl;
}

int main(int argc, char*argv[])
{
        // create FFT class implementation
        // each implementation has its own plan map.
        FFT<double> fft;
        // create some data
        size_t N = 1024 * 1024;
        std::vector<complex<double>,fftalloc<complex<double> > > data(N);
        std::vector<complex<double>,fftalloc<complex<double> > > dataFourier(N);
        std::vector<complex<double>,fftalloc<complex<double> > > dataCalc(N);
        std::cout << "1) --- original data ---" << endl;
        for(std::vector<complex<double> >::size_type i = 0; i != data.size(); i++) {
            data[i] = polar(sin(double(i)/N * M_PI * 2), 0.0);
//            std::cout << i << data[i] << endl;
        }

        // fft
        std::cout << "2) --- fft data ---" << endl;
        fft.fwd(dataFourier, data);
        std::cout << "speed / ms: " << fft.speed()/1000 << endl;
//        for(std::vector<complex<double> >::size_type i = 0; i != dataFourier.size(); i++) {
//            std::cout << i << dataFourier[i] << endl;
//        }

        // ifft
        std::cout << "3) --- ifft data ---" << endl;
        fft.inv(dataCalc, dataFourier);
        std::cout << "speed / ms: " << fft.speed()/1000 << endl;
//        for(std::vector<complex<double> >::size_type i = 0; i != dataCalc.size(); i++) {
//            std::cout << i << dataCalc[i] << endl;
//        }
//        std::cout << "4) --- comparison data ---" << endl;
//        for(std::vector<complex<double> >::size_type i = 0; i != dataCalc.size(); i++) {
//            std::cout << i << data[i] - dataCalc[i] << endl;
//        }

    std::getchar();
    return 0;
}

void TMSI::run()
{
    while(m_bIsRunning)
    {
        //std::cout<<"TMSI::run(s)"<<std::endl;

        //pop matrix only if the producer thread is running
        if(m_pTMSIProducer->isRunning())
        {
            MatrixXf matValue = m_pRawMatrixBuffer_In->pop();

            // Set Beep trigger (if activated)
            if(m_bBeepTrigger && m_qTimerTrigger.elapsed() >= m_iTriggerInterval)
            {
                QFuture<void> future = QtConcurrent::run(Beep, 450, 700);
                //Set trigger in received data samples - just for one sample, so that this event is easy to detect
                matValue(136, m_iSamplesPerBlock-1) = 252;
                m_qTimerTrigger.restart();

                Q_UNUSED(future);
            }

            // Set keyboard trigger (if activated and !=0)
            if(m_bUseKeyboardTrigger && m_iTriggerType!=0)
                matValue(136, m_iSamplesPerBlock-1) = m_iTriggerType;

            //Write raw data to fif file
            if(m_bWriteToFile)
                m_pOutfid->write_raw_buffer(matValue.cast<double>(), m_cals);

            // TODO: Use preprocessing if wanted by the user
            if(m_bUseFiltering)
            {
                MatrixXf temp = matValue;

                matValue = matValue - m_matOldMatrix;
                m_matOldMatrix = temp;

                //    //Check filter class - will be removed in the future - testing purpose only!
                //    FilterTools* filterObject = new FilterTools();

                //    //kaiser window testing
                //    qint32 numberCoeff = 51;
                //    QVector<float> impulseResponse(numberCoeff);
                //    filterObject->createDynamicFilter(QString('LP'), numberCoeff, (float)0.3, impulseResponse);

                //    ofstream outputFileStream("mne_x_plugins/resources/tmsi/filterToolsTest.txt", ios::out);

                //    outputFileStream << "impulseResponse:\n";
                //    for(int i=0; i<impulseResponse.size(); i++)
                //        outputFileStream << impulseResponse[i] << " ";
                //    outputFileStream << endl;

                //    //convolution testing
                //    QVector<float> in (12, 2);
                //    QVector<float> kernel (4, 2);

                //    QVector<float> out = filterObject->convolve(in, kernel);

                //    outputFileStream << "convolution result:\n";
                //    for(int i=0; i<out.size(); i++)
                //        outputFileStream << out[i] << " ";
                //    outputFileStream << endl;
            }

            // TODO: Perform a fft if wanted by the user
            if(m_bUseFFT)
            {
                QElapsedTimer timer;
                timer.start();

                FFT<float> fft;
                Matrix<complex<float>, 138, 16> freq;

                for(qint32 i = 0; i < matValue.rows(); ++i)
                    fft.fwd(freq.row(i), matValue.row(i));

//                cout<<"FFT postprocessing done in "<<timer.nsecsElapsed()<<" nanosec"<<endl;
//                cout<<"matValue before FFT:"<<endl<<matValue<<endl;
//                cout<<"freq after FFT:"<<endl<<freq<<endl;
//                matValue = freq.cwiseAbs();
//                cout<<"matValue after FFT:"<<endl<<matValue<<endl;
            }

            //Change values of the trigger channel for better plotting - this change is not saved in the produced fif file
            if(m_iNumberOfChannels>137)
            {
                for(int i = 0; i<matValue.row(137).cols(); i++)
                {
                    // Left keyboard or capacitive
                    if(matValue.row(136)[i] == 254)
                        matValue.row(136)[i] = 4000;

                    // Right keyboard
                    if(matValue.row(136)[i] == 253)
                        matValue.row(136)[i] = 8000;

                    // Beep
                    if(matValue.row(136)[i] == 252)
                        matValue.row(136)[i] = 2000;
//.........这里部分代码省略.........

void UnbiasedSquaredPhaseLagIndex::compute(ConnectivitySettings::IntermediateTrialData& inputData,
                                           QVector<QPair<int,MatrixXcd> >& vecPairCsdSum,
                                           QVector<QPair<int,MatrixXd> >& vecPairCsdImagSignSum,
                                           QMutex& mutex,
                                           int iNRows,
                                           int iNFreqs,
                                           int iNfft,
                                           const QPair<MatrixXd, VectorXd>& tapers)
{
    if(inputData.vecPairCsdImagSign.size() == iNRows) {
        //qDebug() << "UnbiasedSquaredPhaseLagIndex::compute - vecPairCsdImagSign was already computed for this trial.";
        return;
    }

    inputData.vecPairCsdImagSign.clear();

    int i,j;

    // Calculate tapered spectra if not available already
    // This code was copied and changed modified Utils/Spectra since we do not want to call the function due to time loss.
    if(inputData.vecTapSpectra.size() != iNRows) {
        inputData.vecTapSpectra.clear();

        RowVectorXd vecInputFFT, rowData;
        RowVectorXcd vecTmpFreq;

        MatrixXcd matTapSpectrum(tapers.first.rows(), iNFreqs);

        QVector<Eigen::MatrixXcd> vecTapSpectra;

        FFT<double> fft;
        fft.SetFlag(fft.HalfSpectrum);

        for (i = 0; i < iNRows; ++i) {
            // Substract mean
            rowData.array() = inputData.matData.row(i).array() - inputData.matData.row(i).mean();

            // Calculate tapered spectra if not available already
            for(j = 0; j < tapers.first.rows(); j++) {
                vecInputFFT = rowData.cwiseProduct(tapers.first.row(j));
                // FFT for freq domain returning the half spectrum and multiply taper weights
                fft.fwd(vecTmpFreq, vecInputFFT, iNfft);
                matTapSpectrum.row(j) = vecTmpFreq * tapers.second(j);
            }

            inputData.vecTapSpectra.append(matTapSpectrum);
        }
    }

    // Compute CSD
    if(inputData.vecPairCsd.isEmpty()) {
        double denomCSD = sqrt(tapers.second.cwiseAbs2().sum()) * sqrt(tapers.second.cwiseAbs2().sum()) / 2.0;

        bool bNfftEven = false;
        if (iNfft % 2 == 0){
            bNfftEven = true;
        }

        MatrixXcd matCsd = MatrixXcd(iNRows, iNFreqs);

        for (i = 0; i < iNRows; ++i) {
            for (j = i; j < iNRows; ++j) {
                // Compute CSD (average over tapers if necessary)
                matCsd.row(j) = inputData.vecTapSpectra.at(i).cwiseProduct(inputData.vecTapSpectra.at(j).conjugate()).colwise().sum() / denomCSD;

                // Divide first and last element by 2 due to half spectrum
                matCsd.row(j)(0) /= 2.0;
                if(bNfftEven) {
                    matCsd.row(j).tail(1) /= 2.0;
                }
            }

            inputData.vecPairCsd.append(QPair<int,MatrixXcd>(i,matCsd));
            inputData.vecPairCsdImagSign.append(QPair<int,MatrixXd>(i,matCsd.imag().cwiseSign()));
        }

        mutex.lock();

        if(vecPairCsdSum.isEmpty()) {
            vecPairCsdSum = inputData.vecPairCsd;
            vecPairCsdImagSignSum = inputData.vecPairCsdImagSign;
        } else {
            for (int j = 0; j < vecPairCsdSum.size(); ++j) {
                vecPairCsdSum[j].second += inputData.vecPairCsd.at(j).second;
                vecPairCsdImagSignSum[j].second += inputData.vecPairCsdImagSign.at(j).second;
            }
        }

        mutex.unlock();
    } else {
        if(inputData.vecPairCsdImagSign.isEmpty()) {
            for (i = 0; i < inputData.vecPairCsd.size(); ++i) {
                inputData.vecPairCsdImagSign.append(QPair<int,MatrixXd>(i,inputData.vecPairCsd.at(i).second.imag().cwiseSign()));
            }

            mutex.lock();

            if(vecPairCsdImagSignSum.isEmpty()) {
                vecPairCsdImagSignSum = inputData.vecPairCsdImagSign;
            } else {
//.........这里部分代码省略.........

boolean CAlgorithmHilbertTransform::process(void)
{

	uint32 l_ui32ChannelCount = ip_pMatrix->getDimensionSize(0);
	uint32 l_ui32SamplesPerChannel = ip_pMatrix->getDimensionSize(1);

	IMatrix* l_pInputMatrix = ip_pMatrix;
	IMatrix* l_pOutputEnvelopeMatrix = op_pEnvelopeMatrix;
	IMatrix* l_pOutputPhaseMatrix = op_pPhaseMatrix;

	FFT< double, internal::kissfft_impl<double > > fft;

	if(this->isInputTriggerActive(OVP_Algorithm_HilbertTransform_InputTriggerId_Process))
	{
		//Computing Hilbert transform for all channels
		for(uint32 channel=0; channel<l_ui32ChannelCount; channel++)
		{
			//Initialization of buffer vectors
			m_vecXcdSignalBuffer = RowVectorXcd::Zero(l_ui32SamplesPerChannel);
			m_vecXcdSignalFourier = RowVectorXcd::Zero(l_ui32SamplesPerChannel);

			//Initialization of vector h used to compute analytic signal
			m_vecXdHilbert.resize(l_ui32SamplesPerChannel);
			m_vecXdHilbert(0) = 1.0;
			if(l_ui32SamplesPerChannel%2 == 0)
			{
				m_vecXdHilbert(l_ui32SamplesPerChannel/2) = 1.0;
				for(uint32 i=1; i<l_ui32SamplesPerChannel/2; i++)
				{
					m_vecXdHilbert(i) = 2.0;
				}
				for(uint32 i=(l_ui32SamplesPerChannel/2)+1; i<l_ui32SamplesPerChannel; i++)
				{
					m_vecXdHilbert(i) = 0.0;
				}
			}
			else
			{
				m_vecXdHilbert((l_ui32SamplesPerChannel/2)+1) = 1.0;
				for(uint32 i=1; i<(l_ui32SamplesPerChannel/2)+1; i++)
				{
					m_vecXdHilbert(i) = 2.0;
				}
					for(uint32 i=(l_ui32SamplesPerChannel/2)+2; i<l_ui32SamplesPerChannel; i++)
				{
					m_vecXdHilbert(i) = 0.0;
				}
			}

			//Copy input signal chunk on buffer
			for(uint32 samples=0; samples<l_ui32SamplesPerChannel;samples++)
			{
				m_vecXcdSignalBuffer(samples).real(l_pInputMatrix->getBuffer()[samples + channel * (l_ui32SamplesPerChannel)]);
				m_vecXcdSignalBuffer(samples).imag(0.0);
			}

			//Fast Fourier Transform of input signal
			fft.fwd(m_vecXcdSignalFourier, m_vecXcdSignalBuffer);

			//Apply Hilbert transform by element-wise multiplying fft vector by h
			m_vecXcdSignalFourier = m_vecXcdSignalFourier * m_vecXdHilbert;

			//Inverse Fast Fourier transform
			fft.inv(m_vecXcdSignalBuffer, m_vecXcdSignalFourier); //m_vecXcdSignalBuffer is now the analytical signal of the initial input signal

			//Compute envelope and phase and pass it to the corresponding output
			for(uint32 samples=0; samples<l_ui32SamplesPerChannel;samples++)
			{
				l_pOutputEnvelopeMatrix->getBuffer()[samples + channel*samples] = abs(m_vecXcdSignalBuffer(samples));
				l_pOutputPhaseMatrix->getBuffer()[samples + channel*samples] = arg(m_vecXcdSignalBuffer(samples));
			}

		}
	}
	return true;
}

boolean CAlgorithmHilbertTransform::process(void)
{

	uint32 l_ui32ChannelCount = ip_pMatrix->getDimensionSize(0);
	uint32 l_ui32SamplesPerChannel = ip_pMatrix->getDimensionSize(1);


	IMatrix* l_pInputMatrix = ip_pMatrix;
	IMatrix* l_pOutputHilbertMatrix = op_pHilbertMatrix;
	IMatrix* l_pOutputEnvelopeMatrix = op_pEnvelopeMatrix;
	IMatrix* l_pOutputPhaseMatrix = op_pPhaseMatrix;

	FFT< double, internal::kissfft_impl<double > > fft; //create instance of fft transform

	if(this->isInputTriggerActive(OVP_Algorithm_HilbertTransform_InputTriggerId_Initialize)) //Check if the input is correct
	{
		if( l_pInputMatrix->getDimensionCount() != 2)
		{
			this->getLogManager() << LogLevel_Error << "The input matrix must have 2 dimensions, here the dimension is ";
			std::cout<<l_pInputMatrix->getDimensionCount()<<std::endl;
			return false;
		}

		//Setting size of outputs

		l_pOutputHilbertMatrix->setDimensionCount(2);
		l_pOutputHilbertMatrix->setDimensionSize(0,l_ui32ChannelCount);
		l_pOutputHilbertMatrix->setDimensionSize(1,l_ui32SamplesPerChannel);

		l_pOutputEnvelopeMatrix->setDimensionCount(2);
		l_pOutputEnvelopeMatrix->setDimensionSize(0,l_ui32ChannelCount);
		l_pOutputEnvelopeMatrix->setDimensionSize(1,l_ui32SamplesPerChannel);

		l_pOutputPhaseMatrix->setDimensionCount(2);
		l_pOutputPhaseMatrix->setDimensionSize(0,l_ui32ChannelCount);
		l_pOutputPhaseMatrix->setDimensionSize(1,l_ui32SamplesPerChannel);

	}

	if(this->isInputTriggerActive(OVP_Algorithm_HilbertTransform_InputTriggerId_Process))
	{

		//Computing Hilbert transform for all channels
		for(uint32 channel=0; channel<l_ui32ChannelCount; channel++)
		{
			//Initialization of buffer vectors
			m_vecXcdSignalBuffer = VectorXcd::Zero(l_ui32SamplesPerChannel);
			m_vecXcdSignalFourier = VectorXcd::Zero(l_ui32SamplesPerChannel);

			//Initialization of vector h used to compute analytic signal
			m_vecXdHilbert.resize(l_ui32SamplesPerChannel);
			m_vecXdHilbert(0) = 1.0;

			if(l_ui32SamplesPerChannel%2 == 0)
			{
				m_vecXdHilbert(l_ui32SamplesPerChannel/2) = 1.0;
				for(uint32 i=1; i<l_ui32SamplesPerChannel/2; i++)
				{
					m_vecXdHilbert(i) = 2.0;
				}
				for(uint32 i=(l_ui32SamplesPerChannel/2)+1; i<l_ui32SamplesPerChannel; i++)
				{
					m_vecXdHilbert(i) = 0.0;
				}
			}
			else
			{
				m_vecXdHilbert((l_ui32SamplesPerChannel+1)/2) = 1.0;
				for(uint32 i=1; i<(l_ui32SamplesPerChannel+1); i++)
				{
					m_vecXdHilbert(i) = 2.0;
				}
				for(uint32 i=(l_ui32SamplesPerChannel+1)/2+1; i<l_ui32SamplesPerChannel; i++)
				{
					m_vecXdHilbert(i) = 0.0;
				}
			}

			//Copy input signal chunk on buffer
			for(uint32 samples=0; samples<l_ui32SamplesPerChannel;samples++)
			{
				m_vecXcdSignalBuffer(samples).real(l_pInputMatrix->getBuffer()[samples + channel * (l_ui32SamplesPerChannel)]);
				m_vecXcdSignalBuffer(samples).imag(0.0);
			}

			//Fast Fourier Transform of input signal
			fft.fwd(m_vecXcdSignalFourier, m_vecXcdSignalBuffer);

			//Apply Hilbert transform by element-wise multiplying fft vector by h
			for(uint32 samples=0; samples<l_ui32SamplesPerChannel;samples++)
			{
				m_vecXcdSignalFourier(samples) = m_vecXcdSignalFourier(samples)*m_vecXdHilbert(samples);
			}

			//Inverse Fast Fourier transform
			fft.inv(m_vecXcdSignalBuffer, m_vecXcdSignalFourier); // m_vecXcdSignalBuffer is now the analytical signal of the initial input signal

			//Compute envelope and phase and pass it to the corresponding output
			for(uint32 samples=0; samples<l_ui32SamplesPerChannel;samples++)
			{
//.........这里部分代码省略.........