Mercurial > hg > qm-dsp
comparison dsp/rateconversion/Resampler.cpp @ 366:767947956fc1
More resampler fixes (particularly to latency calculation) and tests
author | Chris Cannam <c.cannam@qmul.ac.uk> |
---|---|
date | Mon, 14 Oct 2013 16:15:32 +0100 |
parents | b73dad5e6201 |
children | f8fc21365a8c |
comparison
equal
deleted
inserted
replaced
365:b73dad5e6201 | 366:767947956fc1 |
---|---|
8 | 8 |
9 #include <iostream> | 9 #include <iostream> |
10 #include <vector> | 10 #include <vector> |
11 | 11 |
12 using std::vector; | 12 using std::vector; |
13 | |
14 //#define DEBUG_RESAMPLER 1 | |
13 | 15 |
14 Resampler::Resampler(int sourceRate, int targetRate) : | 16 Resampler::Resampler(int sourceRate, int targetRate) : |
15 m_sourceRate(sourceRate), | 17 m_sourceRate(sourceRate), |
16 m_targetRate(targetRate) | 18 m_targetRate(targetRate) |
17 { | 19 { |
50 kw.cut(filter); | 52 kw.cut(filter); |
51 | 53 |
52 int inputSpacing = m_targetRate / m_gcd; | 54 int inputSpacing = m_targetRate / m_gcd; |
53 int outputSpacing = m_sourceRate / m_gcd; | 55 int outputSpacing = m_sourceRate / m_gcd; |
54 | 56 |
55 m_latency = int(ceil((m_filterLength / 2.0) / outputSpacing)); | 57 #ifdef DEBUG_RESAMPLER |
56 | 58 std::cerr << "resample " << m_sourceRate << " -> " << m_targetRate |
57 int bufferLength = 0; | 59 << ": inputSpacing " << inputSpacing << ", outputSpacing " |
60 << outputSpacing << ": filter length " << m_filterLength | |
61 << std::endl; | |
62 #endif | |
58 | 63 |
59 m_phaseData = new Phase[inputSpacing]; | 64 m_phaseData = new Phase[inputSpacing]; |
60 | 65 |
61 for (int phase = 0; phase < inputSpacing; ++phase) { | 66 for (int phase = 0; phase < inputSpacing; ++phase) { |
62 | 67 |
64 | 69 |
65 p.nextPhase = phase - outputSpacing; | 70 p.nextPhase = phase - outputSpacing; |
66 while (p.nextPhase < 0) p.nextPhase += inputSpacing; | 71 while (p.nextPhase < 0) p.nextPhase += inputSpacing; |
67 p.nextPhase %= inputSpacing; | 72 p.nextPhase %= inputSpacing; |
68 | 73 |
69 p.drop = int(ceil(std::max(0, outputSpacing - phase) / inputSpacing)); | 74 p.drop = int(ceil(std::max(0.0, double(outputSpacing - phase)) |
70 p.take = int((outputSpacing + | 75 / inputSpacing)); |
71 ((m_filterLength - 1 - phase) % inputSpacing)) | 76 |
72 / outputSpacing); | 77 int filtZipLength = int(ceil(double(m_filterLength - phase) |
73 | 78 / inputSpacing)); |
74 int filtZipLength = int(ceil((m_filterLength - phase) / inputSpacing)); | |
75 if (filtZipLength > bufferLength) { | |
76 bufferLength = filtZipLength; | |
77 } | |
78 | |
79 for (int i = 0; i < filtZipLength; ++i) { | 79 for (int i = 0; i < filtZipLength; ++i) { |
80 p.filter.push_back(filter[i * inputSpacing + phase]); | 80 p.filter.push_back(filter[i * inputSpacing + phase]); |
81 } | 81 } |
82 | 82 |
83 m_phaseData[phase] = p; | 83 m_phaseData[phase] = p; |
84 } | 84 } |
85 | |
86 #ifdef DEBUG_RESAMPLER | |
87 for (int phase = 0; phase < inputSpacing; ++phase) { | |
88 std::cerr << "filter for phase " << phase << " of " << inputSpacing << " (with length " << m_phaseData[phase].filter.size() << "):"; | |
89 for (int i = 0; i < m_phaseData[phase].filter.size(); ++i) { | |
90 if (i % 4 == 0) { | |
91 std::cerr << std::endl << i << ": "; | |
92 } | |
93 float v = m_phaseData[phase].filter[i]; | |
94 if (v == 1) { | |
95 std::cerr << " *** " << v << " *** "; | |
96 } else { | |
97 std::cerr << v << " "; | |
98 } | |
99 } | |
100 std::cerr << std::endl; | |
101 } | |
102 #endif | |
85 | 103 |
86 delete[] filter; | 104 delete[] filter; |
87 | 105 |
88 // The May implementation of this uses a pull model -- we ask the | 106 // The May implementation of this uses a pull model -- we ask the |
89 // resampler for a certain number of output samples, and it asks | 107 // resampler for a certain number of output samples, and it asks |
100 // either return a very variable number of samples (none at all | 118 // either return a very variable number of samples (none at all |
101 // until the filter fills, then half the filter length at once) or | 119 // until the filter fills, then half the filter length at once) or |
102 // else have a lengthy declared latency on the output. We do the | 120 // else have a lengthy declared latency on the output. We do the |
103 // latter. (What do other implementations do?) | 121 // latter. (What do other implementations do?) |
104 | 122 |
105 m_phase = m_filterLength % inputSpacing; | 123 m_phase = (m_filterLength/2) % inputSpacing; |
106 m_buffer = vector<double>(bufferLength, 0); | 124 |
125 m_buffer = vector<double>(m_phaseData[0].filter.size(), 0); | |
126 | |
127 m_latency = | |
128 ((m_buffer.size() * inputSpacing) - (m_filterLength/2)) / outputSpacing | |
129 + m_phase; | |
130 | |
131 #ifdef DEBUG_RESAMPLER | |
132 std::cerr << "initial phase " << m_phase << " (as " << (m_filterLength/2) << " % " << inputSpacing << ")" | |
133 << ", latency " << m_latency << std::endl; | |
134 #endif | |
107 } | 135 } |
108 | 136 |
109 double | 137 double |
110 Resampler::reconstructOne(const double *src) | 138 Resampler::reconstructOne() |
111 { | 139 { |
112 Phase &pd = m_phaseData[m_phase]; | 140 Phase &pd = m_phaseData[m_phase]; |
113 double *filt = pd.filter.data(); | 141 double *filt = pd.filter.data(); |
142 double v = 0.0; | |
114 int n = pd.filter.size(); | 143 int n = pd.filter.size(); |
115 double v = 0.0; | |
116 for (int i = 0; i < n; ++i) { | 144 for (int i = 0; i < n; ++i) { |
117 v += m_buffer[i] * filt[i]; | 145 v += m_buffer[i] * filt[i]; |
118 } | 146 } |
119 m_buffer = vector<double>(m_buffer.begin() + pd.drop, m_buffer.end()); | 147 m_buffer = vector<double>(m_buffer.begin() + pd.drop, m_buffer.end()); |
120 for (int i = 0; i < pd.take; ++i) { | 148 m_phase = pd.nextPhase; |
149 return v; | |
150 } | |
151 | |
152 int | |
153 Resampler::process(const double *src, double *dst, int n) | |
154 { | |
155 for (int i = 0; i < n; ++i) { | |
121 m_buffer.push_back(src[i]); | 156 m_buffer.push_back(src[i]); |
122 } | 157 } |
123 return v; | 158 |
124 } | 159 int maxout = int(ceil(double(n) * m_targetRate / m_sourceRate)); |
125 | 160 int outidx = 0; |
126 int | 161 |
127 Resampler::process(const double *src, double *dst, int remaining) | 162 #ifdef DEBUG_RESAMPLER |
128 { | 163 std::cerr << "process: buf siz " << m_buffer.size() << " filt siz for phase " << m_phase << " " << m_phaseData[m_phase].filter.size() << std::endl; |
129 int m = 0; | 164 #endif |
130 int offset = 0; | 165 |
131 | 166 while (outidx < maxout && |
132 while (remaining >= m_phaseData[m_phase].take) { | 167 m_buffer.size() >= m_phaseData[m_phase].filter.size()) { |
133 // std::cerr << "remaining = " << remaining << ", m = " << m << ", take = " << m_phaseData[m_phase].take << std::endl; | 168 dst[outidx] = reconstructOne(); |
134 int advance = m_phaseData[m_phase].take; | 169 outidx++; |
135 dst[m] = reconstructOne(src + offset); | 170 } |
136 offset += advance; | 171 |
137 remaining -= advance; | 172 return outidx; |
138 m_phase = m_phaseData[m_phase].nextPhase; | 173 } |
139 // std::cerr << "remaining -> " << remaining << ", new phase has advance " << m_phaseData[m_phase].take << std::endl; | 174 |
140 ++m; | |
141 } | |
142 | |
143 // if (remaining > 0) { | |
144 // std::cerr << "have " << remaining << " spare, pushing to buffer" << std::endl; | |
145 // } | |
146 | |
147 for (int i = 0; i < remaining; ++i) { | |
148 m_buffer.push_back(src[offset + i]); | |
149 } | |
150 | |
151 return m; | |
152 } | |
153 | |
154 std::vector<double> | 175 std::vector<double> |
155 Resampler::resample(int sourceRate, int targetRate, const double *data, int n) | 176 Resampler::resample(int sourceRate, int targetRate, const double *data, int n) |
156 { | 177 { |
157 Resampler r(sourceRate, targetRate); | 178 Resampler r(sourceRate, targetRate); |
158 | 179 |
159 int latency = r.getLatency(); | 180 int latency = r.getLatency(); |
160 | 181 |
161 int m = int(ceil((n * targetRate) / sourceRate)); | 182 int m = int(ceil(double(n * targetRate) / sourceRate)); |
162 int m1 = m + latency; | 183 int m1 = m + latency; |
163 int n1 = int((m1 * sourceRate) / targetRate); | 184 int n1 = int(double(m1 * sourceRate) / targetRate); |
164 | 185 |
165 vector<double> pad(n1 - n, 0.0); | 186 vector<double> pad(n1 - n, 0.0); |
166 vector<double> out(m1, 0.0); | 187 vector<double> out(m1, 0.0); |
167 | 188 |
168 int got = r.process(data, out.data(), n); | 189 int got = r.process(data, out.data(), n); |
169 got += r.process(pad.data(), out.data() + got, pad.size()); | 190 got += r.process(pad.data(), out.data() + got, pad.size()); |
170 | 191 |
192 #ifdef DEBUG_RESAMPLER | |
193 std::cerr << "resample: " << n << " in, " << got << " out" << std::endl; | |
194 for (int i = 0; i < got; ++i) { | |
195 if (i % 5 == 0) std::cout << std::endl << i << "... "; | |
196 std::cout << (float) out[i] << " "; | |
197 } | |
198 std::cout << std::endl; | |
199 #endif | |
200 | |
171 return vector<double>(out.begin() + latency, out.begin() + got); | 201 return vector<double>(out.begin() + latency, out.begin() + got); |
172 } | 202 } |
173 | 203 |