FFTRealPassInverse.hpp (6231B)
1 /***************************************************************************** 2 3 FFTRealPassInverse.hpp 4 By Laurent de Soras 5 6 --- Legal stuff --- 7 8 This program is free software. It comes without any warranty, to 9 the extent permitted by applicable law. You can redistribute it 10 and/or modify it under the terms of the Do What The Fuck You Want 11 To Public License, Version 2, as published by Sam Hocevar. See 12 http://sam.zoy.org/wtfpl/COPYING for more details. 13 14 *Tab=3***********************************************************************/ 15 16 17 18 #if defined (ffft_FFTRealPassInverse_CURRENT_CODEHEADER) 19 #error Recursive inclusion of FFTRealPassInverse code header. 20 #endif 21 #define ffft_FFTRealPassInverse_CURRENT_CODEHEADER 22 23 #if ! defined (ffft_FFTRealPassInverse_CODEHEADER_INCLUDED) 24 #define ffft_FFTRealPassInverse_CODEHEADER_INCLUDED 25 26 27 28 /*\\\ INCLUDE FILES \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ 29 30 #include "ffft/FFTRealUseTrigo.h" 31 32 33 34 namespace ffft 35 { 36 37 38 39 /*\\\ PUBLIC \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ 40 41 42 43 template <int PASS> 44 void FFTRealPassInverse <PASS>::process (long len, DataType dest_ptr [], DataType src_ptr [], const DataType f_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 45 { 46 process_internal ( 47 len, 48 dest_ptr, 49 f_ptr, 50 cos_ptr, 51 cos_len, 52 br_ptr, 53 osc_list 54 ); 55 FFTRealPassInverse <PASS - 1>::process_rec ( 56 len, 57 src_ptr, 58 dest_ptr, 59 cos_ptr, 60 cos_len, 61 br_ptr, 62 osc_list 63 ); 64 } 65 66 67 68 template <int PASS> 69 void FFTRealPassInverse <PASS>::process_rec (long len, DataType dest_ptr [], DataType src_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 70 { 71 process_internal ( 72 len, 73 dest_ptr, 74 src_ptr, 75 cos_ptr, 76 cos_len, 77 br_ptr, 78 osc_list 79 ); 80 FFTRealPassInverse <PASS - 1>::process_rec ( 81 len, 82 src_ptr, 83 dest_ptr, 84 cos_ptr, 85 cos_len, 86 br_ptr, 87 osc_list 88 ); 89 } 90 91 template <> 92 inline void FFTRealPassInverse <0>::process_rec (long len, DataType dest_ptr [], DataType src_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 93 { 94 // Stops recursion 95 } 96 97 98 99 template <int PASS> 100 void FFTRealPassInverse <PASS>::process_internal (long len, DataType dest_ptr [], const DataType src_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 101 { 102 const long dist = 1L << (PASS - 1); 103 const long c1_r = 0; 104 const long c1_i = dist; 105 const long c2_r = dist * 2; 106 const long c2_i = dist * 3; 107 const long cend = dist * 4; 108 const long table_step = cos_len >> (PASS - 1); 109 110 enum { TRIGO_OSC = PASS - FFTRealFixLenParam::TRIGO_BD_LIMIT }; 111 enum { TRIGO_DIRECT = (TRIGO_OSC >= 0) ? 1 : 0 }; 112 113 long coef_index = 0; 114 do 115 { 116 const DataType * const sf = src_ptr + coef_index; 117 DataType * const df = dest_ptr + coef_index; 118 119 // Extreme coefficients are always real 120 df [c1_r] = sf [c1_r] + sf [c2_r]; 121 df [c2_r] = sf [c1_r] - sf [c2_r]; 122 df [c1_i] = sf [c1_i] * 2; 123 df [c2_i] = sf [c2_i] * 2; 124 125 FFTRealUseTrigo <TRIGO_DIRECT>::prepare (osc_list [TRIGO_OSC]); 126 127 // Others are conjugate complex numbers 128 for (long i = 1; i < dist; ++ i) 129 { 130 df [c1_r + i] = sf [c1_r + i] + sf [c2_r - i]; 131 df [c1_i + i] = sf [c2_r + i] - sf [cend - i]; 132 133 DataType c; 134 DataType s; 135 FFTRealUseTrigo <TRIGO_DIRECT>::iterate ( 136 osc_list [TRIGO_OSC], 137 c, 138 s, 139 cos_ptr, 140 i * table_step, 141 (dist - i) * table_step 142 ); 143 144 const DataType vr = sf [c1_r + i] - sf [c2_r - i]; 145 const DataType vi = sf [c2_r + i] + sf [cend - i]; 146 147 df [c2_r + i] = vr * c + vi * s; 148 df [c2_i + i] = vi * c - vr * s; 149 } 150 151 coef_index += cend; 152 } 153 while (coef_index < len); 154 } 155 156 template <> 157 inline void FFTRealPassInverse <2>::process_internal (long len, DataType dest_ptr [], const DataType src_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 158 { 159 // Antepenultimate pass 160 const DataType sqrt2_2 = DataType (SQRT2 * 0.5); 161 162 long coef_index = 0; 163 do 164 { 165 dest_ptr [coef_index ] = src_ptr [coef_index] + src_ptr [coef_index + 4]; 166 dest_ptr [coef_index + 4] = src_ptr [coef_index] - src_ptr [coef_index + 4]; 167 dest_ptr [coef_index + 2] = src_ptr [coef_index + 2] * 2; 168 dest_ptr [coef_index + 6] = src_ptr [coef_index + 6] * 2; 169 170 dest_ptr [coef_index + 1] = src_ptr [coef_index + 1] + src_ptr [coef_index + 3]; 171 dest_ptr [coef_index + 3] = src_ptr [coef_index + 5] - src_ptr [coef_index + 7]; 172 173 const DataType vr = src_ptr [coef_index + 1] - src_ptr [coef_index + 3]; 174 const DataType vi = src_ptr [coef_index + 5] + src_ptr [coef_index + 7]; 175 176 dest_ptr [coef_index + 5] = (vr + vi) * sqrt2_2; 177 dest_ptr [coef_index + 7] = (vi - vr) * sqrt2_2; 178 179 coef_index += 8; 180 } 181 while (coef_index < len); 182 } 183 184 template <> 185 inline void FFTRealPassInverse <1>::process_internal (long len, DataType dest_ptr [], const DataType src_ptr [], const DataType cos_ptr [], long cos_len, const long br_ptr [], OscType osc_list []) 186 { 187 // Penultimate and last pass at once 188 const long qlen = len >> 2; 189 190 long coef_index = 0; 191 do 192 { 193 const long ri_0 = br_ptr [coef_index >> 2]; 194 195 const DataType b_0 = src_ptr [coef_index ] + src_ptr [coef_index + 2]; 196 const DataType b_2 = src_ptr [coef_index ] - src_ptr [coef_index + 2]; 197 const DataType b_1 = src_ptr [coef_index + 1] * 2; 198 const DataType b_3 = src_ptr [coef_index + 3] * 2; 199 200 dest_ptr [ri_0 ] = b_0 + b_1; 201 dest_ptr [ri_0 + 2 * qlen] = b_0 - b_1; 202 dest_ptr [ri_0 + 1 * qlen] = b_2 + b_3; 203 dest_ptr [ri_0 + 3 * qlen] = b_2 - b_3; 204 205 coef_index += 4; 206 } 207 while (coef_index < len); 208 } 209 210 211 212 /*\\\ PROTECTED \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ 213 214 215 216 /*\\\ PRIVATE \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ 217 218 219 220 } // namespace ffft 221 222 223 224 #endif // ffft_FFTRealPassInverse_CODEHEADER_INCLUDED 225 226 #undef ffft_FFTRealPassInverse_CURRENT_CODEHEADER 227 228 229 230 /*\\\ EOF \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/