@@ -26,6 +26,14 @@ namespace yup
2626
2727// ==============================================================================
2828
29+ /* *
30+ Calculate the magnitude response of a filter.
31+
32+ @param filter The filter to calculate the magnitude response of.
33+ @param buffer The buffer to store the magnitude response in.
34+ @param minFreq The minimum frequency to calculate the response at.
35+ @param maxFreq The maximum frequency to calculate the response at.
36+ */
2937template <typename FloatType, typename FilterType>
3038void calculateFilterMagnitudeResponse (FilterType& filter, Span<Complex<FloatType>> buffer, double minFreq, double maxFreq)
3139{
@@ -36,10 +44,9 @@ void calculateFilterMagnitudeResponse (FilterType& filter, Span<Complex<FloatTyp
3644 const double freq = minFreq * std::pow (maxFreq / minFreq, ratio);
3745
3846 // Get complex response
39- auto response = filter.getComplexResponse (freq);
47+ auto magnitude = filter.getMagnitudeResponse (freq);
4048
4149 // Calculate magnitude in dB
42- double magnitude = std::abs (response);
4350 double magnitudeDb = 20.0 * std::log10 (yup::jmax (magnitude, 1e-12 ));
4451
4552 buffer[i] = { static_cast <FloatType> (freq), static_cast <FloatType> (magnitudeDb) };
@@ -48,6 +55,15 @@ void calculateFilterMagnitudeResponse (FilterType& filter, Span<Complex<FloatTyp
4855
4956// ==============================================================================
5057
58+ /* *
59+ Calculate the phase response of a filter.
60+
61+ @param filter The filter to calculate the phase response of.
62+ @param buffer The buffer to store the phase response in.
63+ @param minFreq The minimum frequency to calculate the response at.
64+ @param maxFreq The maximum frequency to calculate the response at.
65+ */
66+
5167template <typename FloatType, typename FilterType>
5268void calculateFilterPhaseResponse (FilterType& filter, Span<Complex<FloatType>> buffer, double minFreq, double maxFreq)
5369{
@@ -58,10 +74,9 @@ void calculateFilterPhaseResponse (FilterType& filter, Span<Complex<FloatType>>
5874 const double freq = minFreq * std::pow (maxFreq / minFreq, ratio);
5975
6076 // Get complex response
61- auto response = filter.getComplexResponse (freq);
77+ auto phaseRad = filter.getPhaseResponse (freq);
6278
6379 // Calculate phase in degrees
64- double phaseRad = std::arg (response);
6580 double phaseDeg = phaseRad * 180.0 / yup::MathConstants<double >::pi;
6681
6782 buffer[i] = { static_cast <FloatType> (freq), static_cast <FloatType> (phaseDeg) };
@@ -70,6 +85,15 @@ void calculateFilterPhaseResponse (FilterType& filter, Span<Complex<FloatType>>
7085
7186// ==============================================================================
7287
88+ /* *
89+ Calculate the group delay of a filter.
90+
91+ @param filter The filter to calculate the group delay of.
92+ @param buffer The buffer to store the group delay in.
93+ @param minFreq The minimum frequency to calculate the response at.
94+ @param maxFreq The maximum frequency to calculate the response at.
95+ @param sampleRate The sample rate of the filter.
96+ */
7397template <typename FloatType, typename FilterType>
7498void calculateFilterGroupDelay (FilterType& filter, Span<Complex<FloatType>> buffer, double minFreq, double maxFreq, double sampleRate)
7599{
@@ -78,26 +102,23 @@ void calculateFilterGroupDelay (FilterType& filter, Span<Complex<FloatType>> buf
78102 // Logarithmic frequency sweep
79103 const double ratio = static_cast <double > (i) / (buffer.size () - 1 );
80104 const double freq = minFreq * std::pow (maxFreq / minFreq, ratio);
105+ const double deltaFreq = freq * 0.01 ; // Small frequency step
81106
82107 // Calculate group delay (numerical derivative of phase)
83108 double groupDelay = 0.0 ;
84109 if (i > 0 && i < buffer.size () - 1 )
85110 {
86- const double deltaFreq = freq * 0.01 ; // Small frequency step
87- auto responseLow = filter.getComplexResponse (freq - deltaFreq);
88- auto responseHigh = filter.getComplexResponse (freq + deltaFreq);
89-
90- double phaseLow = std::arg (responseLow);
91- double phaseHigh = std::arg (responseHigh);
111+ auto phaseLow = filter.getPhaseResponse (freq - deltaFreq);
112+ auto phaseHigh = filter.getPhaseResponse (freq + deltaFreq);
92113
93114 // Unwrap phase difference
94115 double phaseDiff = phaseHigh - phaseLow;
95116 while (phaseDiff > yup::MathConstants<double >::pi)
96- phaseDiff -= 2.0 * yup::MathConstants<double >::pi ;
117+ phaseDiff -= yup::MathConstants<double >::twoPi ;
97118 while (phaseDiff < -yup::MathConstants<double >::pi)
98- phaseDiff += 2.0 * yup::MathConstants<double >::pi ;
119+ phaseDiff += yup::MathConstants<double >::twoPi ;
99120
100- groupDelay = -phaseDiff / (2.0 * deltaFreq * 2.0 * yup::MathConstants<double >::pi ) * sampleRate;
121+ groupDelay = -phaseDiff / (2.0 * deltaFreq * yup::MathConstants<double >::twoPi ) * sampleRate;
101122 }
102123
103124 buffer[i] = { static_cast <FloatType> (freq), static_cast <FloatType> (groupDelay) };
0 commit comments