X-Git-Url: http://git.kpe.io/?p=ctsim.git;a=blobdiff_plain;f=libctsim%2Fprocsignal.cpp;h=41895ac06a37308f297e000bf8cf9b1ad5b43551;hp=a7933b21eb2a4ffb508c88d5caa35ffc834d3723;hb=9b2bb510160bdb56f04847f5b55ab61dd8a47976;hpb=5c6b29ab4885308cc3381af6e0a68f4804956d2e diff --git a/libctsim/procsignal.cpp b/libctsim/procsignal.cpp index a7933b2..41895ac 100644 --- a/libctsim/procsignal.cpp +++ b/libctsim/procsignal.cpp @@ -9,7 +9,7 @@ ** This is part of the CTSim program ** Copyright (C) 1983-2000 Kevin Rosenberg ** -** $Id: procsignal.cpp,v 1.11 2000/12/29 15:45:06 kevin Exp $ +** $Id: procsignal.cpp,v 1.13 2001/01/02 05:34:57 kevin Exp $ ** ** This program is free software; you can redistribute it and/or modify ** it under the terms of the GNU General Public License (version 2) as @@ -189,7 +189,7 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw pEZPlot->plot (pSGP); } #endif - shuffleNaturalToFourierOrder (adFrequencyFilter, m_nFilterPoints); + Fourier::shuffleNaturalToFourierOrder (adFrequencyFilter, m_nFilterPoints); #ifdef HAVE_SGP if (pEZPlot && m_traceLevel >= Trace::TRACE_PLOT) { pEZPlot->ezset ("title Filter Response: Fourier Order"); @@ -199,7 +199,7 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw pEZPlot->plot (pSGP); } #endif - ProcessSignal::finiteFourierTransform (adFrequencyFilter, m_adFilter, m_nFilterPoints, -1); + ProcessSignal::finiteFourierTransform (adFrequencyFilter, m_adFilter, m_nFilterPoints, FORWARD); delete adFrequencyFilter; #ifdef HAVE_SGP if (pEZPlot && m_traceLevel >= Trace::TRACE_PLOT) { @@ -210,7 +210,7 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw pEZPlot->plot (pSGP); } #endif - shuffleFourierToNaturalOrder (m_adFilter, m_nFilterPoints); + Fourier::shuffleFourierToNaturalOrder (m_adFilter, m_nFilterPoints); #ifdef HAVE_SGP if (pEZPlot && m_traceLevel >= Trace::TRACE_PLOT) { pEZPlot->ezset ("title Inverse Fourier Frequency: Natural Order"); @@ -304,7 +304,7 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw pEZPlot->plot (pSGP); } #endif - shuffleNaturalToFourierOrder (m_adFilter, m_nFilterPoints); + Fourier::shuffleNaturalToFourierOrder (m_adFilter, m_nFilterPoints); #ifdef HAVE_SGP if (pEZPlot && m_traceLevel >= Trace::TRACE_PLOT) { pEZPlot->ezset ("title Filter Filter: Fourier Order"); @@ -370,10 +370,10 @@ ProcessSignal::init (const int idFilter, const int idFilterMethod, double dBandw m_adFilter = new double [m_nFilterPoints]; std::complex* acInverseFilter = new std::complex [m_nFilterPoints]; - finiteFourierTransform (adSpatialFilter, acInverseFilter, m_nFilterPoints, 1); + finiteFourierTransform (adSpatialFilter, acInverseFilter, m_nFilterPoints, BACKWARD); delete adSpatialFilter; for (i = 0; i < m_nFilterPoints; i++) - m_adFilter[i] = std::abs(acInverseFilter[i]) * m_dSignalInc; + m_adFilter[i] = std::abs (acInverseFilter[i]) * m_dSignalInc; delete acInverseFilter; #ifdef HAVE_SGP if (pEZPlot && m_traceLevel >= Trace::TRACE_PLOT) { @@ -553,12 +553,12 @@ ProcessSignal::filterSignal (const float constInput[], double output[]) const for (i = m_nSignalPoints; i < m_nFilterPoints; i++) inputSignal[i] = 0; // zeropad std::complex* fftSignal = new std::complex [m_nFilterPoints]; - finiteFourierTransform (inputSignal, fftSignal, m_nFilterPoints, -1); + finiteFourierTransform (inputSignal, fftSignal, m_nFilterPoints, FORWARD); delete inputSignal; for (i = 0; i < m_nFilterPoints; i++) fftSignal[i] *= m_adFilter[i]; double* inverseFourier = new double [m_nFilterPoints]; - finiteFourierTransform (fftSignal, inverseFourier, m_nFilterPoints, 1); + finiteFourierTransform (fftSignal, inverseFourier, m_nFilterPoints, BACKWARD); delete fftSignal; for (i = 0; i < m_nSignalPoints; i++) output[i] = inverseFourier[i]; @@ -570,12 +570,12 @@ ProcessSignal::filterSignal (const float constInput[], double output[]) const for (i = m_nSignalPoints; i < m_nFilterPoints; i++) inputSignal[i] = 0; // zeropad std::complex* fftSignal = new std::complex [m_nFilterPoints]; - finiteFourierTransform (inputSignal, fftSignal, -1); + finiteFourierTransform (inputSignal, fftSignal, FORWARD); delete inputSignal; for (i = 0; i < m_nFilterPoints; i++) fftSignal[i] *= m_adFilter[i]; double* inverseFourier = new double [m_nFilterPoints]; - finiteFourierTransform (fftSignal, inverseFourier, 1); + finiteFourierTransform (fftSignal, inverseFourier, BACKWARD); delete fftSignal; for (i = 0; i < m_nSignalPoints; i++) output[i] = inverseFourier[i]; @@ -733,8 +733,8 @@ ProcessSignal::finiteFourierTransform (const std::complex input[], std:: std::complex sum (0,0); for (int j = 0; j < n; j++) { double angle = i * j * angleIncrement; - std::complex exponentTerm (cos(angle), sin(angle)); - sum += input[j] * exponentTerm; + std::complex exponentTerm (cos(angle), sin(angle)); + sum += input[j] * exponentTerm; } if (direction < 0) { sum /= n; @@ -855,119 +855,3 @@ ProcessSignal::finiteFourierTransform (const std::complex input[], doubl } } -// Odd Number of Points -// Natural Frequency Order: -(n-1)/2...-1,0,1...(n-1)/2 -// Fourier Frequency Order: 0, 1..(n-1)/2,-(n-1)/2...-1 -// Even Number of Points -// Natural Frequency Order: -n/2...-1,0,1...((n/2)-1) -// Fourier Frequency Order: 0,1...((n/2)-1),-n/2...-1 - -void -ProcessSignal::shuffleNaturalToFourierOrder (double* pdVector, const int n) -{ - double* pdTemp = new double [n]; - int i; - if (n % 2) { // Odd - int iHalfN = (n - 1) / 2; - - pdTemp[0] = pdVector[iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1] = pdVector[i + 1 + iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + iHalfN + 1] = pdVector[i]; - } else { // Even - int iHalfN = n / 2; - pdTemp[0] = pdVector[iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1] = pdVector[i + iHalfN]; - for (i = 0; i < iHalfN - 1; i++) - pdTemp[i + iHalfN + 1] = pdVector[i]; - } - - for (i = 0; i < n; i++) - pdVector[i] = pdTemp[i]; - delete pdTemp; -} - -void -ProcessSignal::shuffleNaturalToFourierOrder (std::complex* pdVector, const int n) -{ - std::complex* pdTemp = new std::complex [n]; - int i; - if (n % 2) { // Odd - int iHalfN = (n - 1) / 2; - - pdTemp[0] = pdVector[iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1] = pdVector[i + 1 + iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + iHalfN + 1] = pdVector[i]; - } else { // Even - int iHalfN = n / 2; - pdTemp[0] = pdVector[iHalfN]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1] = pdVector[i + iHalfN]; - for (i = 0; i < iHalfN - 1; i++) - pdTemp[i + iHalfN + 1] = pdVector[i]; - } - - for (i = 0; i < n; i++) - pdVector[i] = pdTemp[i]; - delete [] pdTemp; -} - - -void -ProcessSignal::shuffleFourierToNaturalOrder (double* pdVector, const int n) -{ - double* pdTemp = new double [n]; - int i; - if (n % 2) { // Odd - int iHalfN = (n - 1) / 2; - - pdTemp[iHalfN] = pdVector[0]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1 + iHalfN] = pdVector[i + 1]; - for (i = 0; i < iHalfN; i++) - pdTemp[i] = pdVector[i + iHalfN + 1]; - } else { // Even - int iHalfN = n / 2; - pdTemp[iHalfN] = pdVector[0]; - for (i = 0; i < iHalfN; i++) - pdTemp[i] = pdVector[i + iHalfN]; - for (i = 0; i < iHalfN - 1; i++) - pdTemp[i + iHalfN + 1] = pdVector[i+1]; - } - - for (i = 0; i < n; i++) - pdVector[i] = pdTemp[i]; - delete pdTemp; -} - -void -ProcessSignal::shuffleFourierToNaturalOrder (std::complex* pdVector, const int n) -{ - std::complex* pdTemp = new std::complex [n]; - int i; - if (n % 2) { // Odd - int iHalfN = (n - 1) / 2; - - pdTemp[iHalfN] = pdVector[0]; - for (i = 0; i < iHalfN; i++) - pdTemp[i + 1 + iHalfN] = pdVector[i + 1]; - for (i = 0; i < iHalfN; i++) - pdTemp[i] = pdVector[i + iHalfN + 1]; - } else { // Even - int iHalfN = n / 2; - pdTemp[iHalfN] = pdVector[0]; - for (i = 0; i < iHalfN; i++) - pdTemp[i] = pdVector[i + iHalfN]; - for (i = 0; i < iHalfN - 1; i++) - pdTemp[i + iHalfN + 1] = pdVector[i+1]; - } - - for (i = 0; i < n; i++) - pdVector[i] = pdTemp[i]; - delete [] pdTemp; -} -