Mercurial > hg > segmenter-vamp-plugin
comparison armadillo-3.900.4/include/armadillo_bits/op_fft_meat.hpp @ 49:1ec0e2823891
Switch to using subrepo copies of qm-dsp, nnls-chroma, vamp-plugin-sdk; update Armadillo version; assume build without external BLAS/LAPACK
author | Chris Cannam |
---|---|
date | Thu, 13 Jun 2013 10:25:24 +0100 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
48:69251e11a913 | 49:1ec0e2823891 |
---|---|
1 // Copyright (C) 2013 Conrad Sanderson | |
2 // Copyright (C) 2013 NICTA (www.nicta.com.au) | |
3 // | |
4 // This Source Code Form is subject to the terms of the Mozilla Public | |
5 // License, v. 2.0. If a copy of the MPL was not distributed with this | |
6 // file, You can obtain one at http://mozilla.org/MPL/2.0/. | |
7 | |
8 | |
9 | |
10 //! \addtogroup op_fft | |
11 //! @{ | |
12 | |
13 | |
14 | |
15 // | |
16 // op_fft_real | |
17 | |
18 | |
19 | |
20 template<typename T1> | |
21 inline | |
22 void | |
23 op_fft_real::apply( Mat< std::complex<typename T1::pod_type> >& out, const mtOp<std::complex<typename T1::pod_type>,T1,op_fft_real>& in ) | |
24 { | |
25 arma_extra_debug_sigprint(); | |
26 | |
27 typedef typename T1::pod_type in_eT; | |
28 typedef typename std::complex<in_eT> out_eT; | |
29 | |
30 const Proxy<T1> P(in.m); | |
31 | |
32 const uword n_rows = P.get_n_rows(); | |
33 const uword n_cols = P.get_n_cols(); | |
34 const uword n_elem = P.get_n_elem(); | |
35 | |
36 const bool is_vec = ( (n_rows == 1) || (n_cols == 1) ); | |
37 | |
38 const uword N_orig = (is_vec) ? n_elem : n_rows; | |
39 const uword N_user = (in.aux_uword_b == 0) ? in.aux_uword_a : N_orig; | |
40 | |
41 fft_engine<out_eT,false> worker(N_user); | |
42 | |
43 // no need to worry about aliasing, as we're going from a real object to complex complex, which by definition cannot alias | |
44 | |
45 if(is_vec) | |
46 { | |
47 (n_cols == 1) ? out.set_size(N_user, 1) : out.set_size(1, N_user); | |
48 | |
49 if( (out.n_elem == 0) || (N_orig == 0) ) | |
50 { | |
51 out.zeros(); | |
52 return; | |
53 } | |
54 | |
55 if( (N_user == 1) && (N_orig >= 1) ) | |
56 { | |
57 out[0] = out_eT( P[0] ); | |
58 return; | |
59 } | |
60 | |
61 podarray<out_eT> data(N_user); | |
62 | |
63 out_eT* data_mem = data.memptr(); | |
64 | |
65 if(N_user > N_orig) { arrayops::inplace_set( &data_mem[N_orig], out_eT(0), (N_user - N_orig) ); } | |
66 | |
67 const uword N = (std::min)(N_user, N_orig); | |
68 | |
69 if(Proxy<T1>::prefer_at_accessor == false) | |
70 { | |
71 typename Proxy<T1>::ea_type X = P.get_ea(); | |
72 | |
73 for(uword i=0; i < N; ++i) { data_mem[i] = out_eT( X[i], in_eT(0) ); } | |
74 } | |
75 else | |
76 { | |
77 if(n_cols == 1) | |
78 { | |
79 for(uword i=0; i < N; ++i) { data_mem[i] = out_eT( P.at(i,0), in_eT(0) ); } | |
80 } | |
81 else | |
82 { | |
83 for(uword i=0; i < N; ++i) { data_mem[i] = out_eT( P.at(0,i), in_eT(0) ); } | |
84 } | |
85 } | |
86 | |
87 worker.run( out.memptr(), data_mem ); | |
88 } | |
89 else | |
90 { | |
91 // process each column seperately | |
92 | |
93 out.set_size(N_user, n_cols); | |
94 | |
95 if( (out.n_elem == 0) || (N_orig == 0) ) | |
96 { | |
97 out.zeros(); | |
98 return; | |
99 } | |
100 | |
101 if( (N_user == 1) && (N_orig >= 1) ) | |
102 { | |
103 for(uword col=0; col < n_cols; ++col) { out.at(0,col) = out_eT( P.at(0,col) ); } | |
104 | |
105 return; | |
106 } | |
107 | |
108 podarray<out_eT> data(N_user); | |
109 | |
110 out_eT* data_mem = data.memptr(); | |
111 | |
112 if(N_user > N_orig) { arrayops::inplace_set( &data_mem[N_orig], out_eT(0), (N_user - N_orig) ); } | |
113 | |
114 const uword N = (std::min)(N_user, N_orig); | |
115 | |
116 for(uword col=0; col < n_cols; ++col) | |
117 { | |
118 for(uword i=0; i < N; ++i) { data_mem[i] = P.at(i, col); } | |
119 | |
120 worker.run( out.colptr(col), data_mem ); | |
121 } | |
122 } | |
123 } | |
124 | |
125 | |
126 | |
127 // | |
128 // op_fft_cx | |
129 | |
130 | |
131 template<typename T1> | |
132 inline | |
133 void | |
134 op_fft_cx::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_fft_cx>& in) | |
135 { | |
136 arma_extra_debug_sigprint(); | |
137 | |
138 typedef typename T1::elem_type eT; | |
139 | |
140 const Proxy<T1> P(in.m); | |
141 | |
142 if(P.is_alias(out) == false) | |
143 { | |
144 op_fft_cx::apply_noalias<T1,false>(out, P, in.aux_uword_a, in.aux_uword_b); | |
145 } | |
146 else | |
147 { | |
148 Mat<eT> tmp; | |
149 | |
150 op_fft_cx::apply_noalias<T1,false>(tmp, P, in.aux_uword_a, in.aux_uword_b); | |
151 | |
152 out.steal_mem(tmp); | |
153 } | |
154 } | |
155 | |
156 | |
157 | |
158 template<typename T1, bool inverse> | |
159 inline | |
160 void | |
161 op_fft_cx::apply_noalias(Mat<typename T1::elem_type>& out, const Proxy<T1>& P, const uword a, const uword b) | |
162 { | |
163 arma_extra_debug_sigprint(); | |
164 | |
165 typedef typename T1::elem_type eT; | |
166 | |
167 const uword n_rows = P.get_n_rows(); | |
168 const uword n_cols = P.get_n_cols(); | |
169 const uword n_elem = P.get_n_elem(); | |
170 | |
171 const bool is_vec = ( (n_rows == 1) || (n_cols == 1) ); | |
172 | |
173 const uword N_orig = (is_vec) ? n_elem : n_rows; | |
174 const uword N_user = (b == 0) ? a : N_orig; | |
175 | |
176 fft_engine<eT,inverse> worker(N_user); | |
177 | |
178 if(is_vec) | |
179 { | |
180 (n_cols == 1) ? out.set_size(N_user, 1) : out.set_size(1, N_user); | |
181 | |
182 if( (out.n_elem == 0) || (N_orig == 0) ) | |
183 { | |
184 out.zeros(); | |
185 return; | |
186 } | |
187 | |
188 if( (N_user == 1) && (N_orig >= 1) ) | |
189 { | |
190 out[0] = P[0]; | |
191 return; | |
192 } | |
193 | |
194 if( (N_user > N_orig) || (is_Mat<typename Proxy<T1>::stored_type>::value == false) ) | |
195 { | |
196 podarray<eT> data(N_user); | |
197 | |
198 eT* data_mem = data.memptr(); | |
199 | |
200 if(N_user > N_orig) { arrayops::inplace_set( &data_mem[N_orig], eT(0), (N_user - N_orig) ); } | |
201 | |
202 op_fft_cx::copy_vec( data_mem, P, (std::min)(N_user, N_orig) ); | |
203 | |
204 worker.run( out.memptr(), data_mem ); | |
205 } | |
206 else | |
207 { | |
208 const unwrap< typename Proxy<T1>::stored_type > tmp(P.Q); | |
209 | |
210 worker.run( out.memptr(), tmp.M.memptr() ); | |
211 } | |
212 } | |
213 else | |
214 { | |
215 // process each column seperately | |
216 | |
217 out.set_size(N_user, n_cols); | |
218 | |
219 if( (out.n_elem == 0) || (N_orig == 0) ) | |
220 { | |
221 out.zeros(); | |
222 return; | |
223 } | |
224 | |
225 if( (N_user == 1) && (N_orig >= 1) ) | |
226 { | |
227 for(uword col=0; col < n_cols; ++col) { out.at(0,col) = P.at(0,col); } | |
228 | |
229 return; | |
230 } | |
231 | |
232 if( (N_user > N_orig) || (is_Mat<typename Proxy<T1>::stored_type>::value == false) ) | |
233 { | |
234 podarray<eT> data(N_user); | |
235 | |
236 eT* data_mem = data.memptr(); | |
237 | |
238 if(N_user > N_orig) { arrayops::inplace_set( &data_mem[N_orig], eT(0), (N_user - N_orig) ); } | |
239 | |
240 const uword N = (std::min)(N_user, N_orig); | |
241 | |
242 for(uword col=0; col < n_cols; ++col) | |
243 { | |
244 for(uword i=0; i < N; ++i) { data_mem[i] = P.at(i, col); } | |
245 | |
246 worker.run( out.colptr(col), data_mem ); | |
247 } | |
248 } | |
249 else | |
250 { | |
251 const unwrap< typename Proxy<T1>::stored_type > tmp(P.Q); | |
252 | |
253 for(uword col=0; col < n_cols; ++col) | |
254 { | |
255 worker.run( out.colptr(col), tmp.M.colptr(col) ); | |
256 } | |
257 } | |
258 } | |
259 | |
260 | |
261 // correct the scaling for the inverse transform | |
262 if(inverse == true) | |
263 { | |
264 typedef typename get_pod_type<eT>::result T; | |
265 | |
266 const T k = T(1) / T(N_user); | |
267 | |
268 eT* out_mem = out.memptr(); | |
269 | |
270 const uword out_n_elem = out.n_elem; | |
271 | |
272 for(uword i=0; i < out_n_elem; ++i) { out_mem[i] *= k; } | |
273 } | |
274 } | |
275 | |
276 | |
277 | |
278 template<typename T1> | |
279 arma_hot | |
280 inline | |
281 void | |
282 op_fft_cx::copy_vec(typename Proxy<T1>::elem_type* dest, const Proxy<T1>& P, const uword N) | |
283 { | |
284 arma_extra_debug_sigprint(); | |
285 | |
286 if(is_Mat< typename Proxy<T1>::stored_type >::value == true) | |
287 { | |
288 op_fft_cx::copy_vec_unwrap(dest, P, N); | |
289 } | |
290 else | |
291 { | |
292 op_fft_cx::copy_vec_proxy(dest, P, N); | |
293 } | |
294 } | |
295 | |
296 | |
297 | |
298 template<typename T1> | |
299 arma_hot | |
300 inline | |
301 void | |
302 op_fft_cx::copy_vec_unwrap(typename Proxy<T1>::elem_type* dest, const Proxy<T1>& P, const uword N) | |
303 { | |
304 arma_extra_debug_sigprint(); | |
305 | |
306 const unwrap< typename Proxy<T1>::stored_type > tmp(P.Q); | |
307 | |
308 arrayops::copy(dest, tmp.M.memptr(), N); | |
309 } | |
310 | |
311 | |
312 | |
313 template<typename T1> | |
314 arma_hot | |
315 inline | |
316 void | |
317 op_fft_cx::copy_vec_proxy(typename Proxy<T1>::elem_type* dest, const Proxy<T1>& P, const uword N) | |
318 { | |
319 arma_extra_debug_sigprint(); | |
320 | |
321 if(Proxy<T1>::prefer_at_accessor == false) | |
322 { | |
323 typename Proxy<T1>::ea_type X = P.get_ea(); | |
324 | |
325 for(uword i=0; i < N; ++i) { dest[i] = X[i]; } | |
326 } | |
327 else | |
328 { | |
329 if(P.get_n_cols() == 1) | |
330 { | |
331 for(uword i=0; i < N; ++i) { dest[i] = P.at(i,0); } | |
332 } | |
333 else | |
334 { | |
335 for(uword i=0; i < N; ++i) { dest[i] = P.at(0,i); } | |
336 } | |
337 } | |
338 } | |
339 | |
340 | |
341 | |
342 // | |
343 // op_ifft_cx | |
344 | |
345 | |
346 template<typename T1> | |
347 inline | |
348 void | |
349 op_ifft_cx::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_ifft_cx>& in) | |
350 { | |
351 arma_extra_debug_sigprint(); | |
352 | |
353 typedef typename T1::elem_type eT; | |
354 | |
355 const Proxy<T1> P(in.m); | |
356 | |
357 if(P.is_alias(out) == false) | |
358 { | |
359 op_fft_cx::apply_noalias<T1,true>(out, P, in.aux_uword_a, in.aux_uword_b); | |
360 } | |
361 else | |
362 { | |
363 Mat<eT> tmp; | |
364 | |
365 op_fft_cx::apply_noalias<T1,true>(tmp, P, in.aux_uword_a, in.aux_uword_b); | |
366 | |
367 out.steal_mem(tmp); | |
368 } | |
369 } | |
370 | |
371 | |
372 | |
373 //! @} |