1 /* 2 This software is part of pffft/pfdsp, a set of simple DSP routines. 3 4 Copyright (c) 2014, Andras Retzler <[email protected]> 5 Copyright (c) 2020 Hayati Ayguen <[email protected]> 6 All rights reserved. 7 8 Redistribution and use in source and binary forms, with or without 9 modification, are permitted provided that the following conditions are met: 10 * Redistributions of source code must retain the above copyright 11 notice, this list of conditions and the following disclaimer. 12 * Redistributions in binary form must reproduce the above copyright 13 notice, this list of conditions and the following disclaimer in the 14 documentation and/or other materials provided with the distribution. 15 * Neither the name of the copyright holder nor the 16 names of its contributors may be used to endorse or promote products 17 derived from this software without specific prior written permission. 18 19 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND 20 ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 21 WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 22 DISCLAIMED. IN NO EVENT SHALL ANDRAS RETZLER BE LIABLE FOR ANY 23 DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 24 (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 25 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND 26 ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 27 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 28 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 29 */ 30 31 #pragma once 32 33 #include <stdio.h> 34 #include <stdint.h> 35 36 #ifdef __cplusplus 37 extern "C" { 38 #endif 39 40 41 /* 42 _____ _ 43 / ____| | | 44 | | ___ _ __ ___ _ __ | | _____ __ 45 | | / _ \| '_ ` _ \| '_ \| |/ _ \ \/ / 46 | |___| (_) | | | | | | |_) | | __/> < 47 \_____\___/|_| |_| |_| .__/|_|\___/_/\_\ 48 | | 49 |_| 50 */ 51 52 typedef struct complexf_s { float i; float q; } complexf; 53 54 // ================================================================================= 55 56 int have_sse_shift_mixer_impl(); 57 58 59 /*********************************************************************/ 60 61 /**************/ 62 /*** ALGO A ***/ 63 /**************/ 64 65 float shift_math_cc(complexf *input, complexf* output, int input_size, float rate, float starting_phase); 66 67 68 /*********************************************************************/ 69 70 /**************/ 71 /*** ALGO B ***/ 72 /**************/ 73 74 typedef struct shift_table_data_s 75 { 76 float* table; 77 int table_size; 78 } shift_table_data_t; 79 80 void shift_table_deinit(shift_table_data_t table_data); 81 shift_table_data_t shift_table_init(int table_size); 82 float shift_table_cc(complexf* input, complexf* output, int input_size, float rate, shift_table_data_t table_data, float starting_phase); 83 84 /*********************************************************************/ 85 86 /**************/ 87 /*** ALGO C ***/ 88 /**************/ 89 90 typedef struct shift_addfast_data_s 91 { 92 float dsin[4]; 93 float dcos[4]; 94 float phase_increment; 95 } shift_addfast_data_t; 96 97 shift_addfast_data_t shift_addfast_init(float rate); 98 float shift_addfast_cc(complexf *input, complexf* output, int input_size, shift_addfast_data_t* d, float starting_phase); 99 float shift_addfast_inp_c(complexf *in_out, int N_cplx, shift_addfast_data_t* d, float starting_phase); 100 101 102 /*********************************************************************/ 103 104 /**************/ 105 /*** ALGO D ***/ 106 /**************/ 107 108 typedef struct shift_unroll_data_s 109 { 110 float* dsin; 111 float* dcos; 112 float phase_increment; 113 int size; 114 } shift_unroll_data_t; 115 116 shift_unroll_data_t shift_unroll_init(float rate, int size); 117 void shift_unroll_deinit(shift_unroll_data_t* d); 118 float shift_unroll_cc(complexf *input, complexf* output, int size, shift_unroll_data_t* d, float starting_phase); 119 float shift_unroll_inp_c(complexf* in_out, int size, shift_unroll_data_t* d, float starting_phase); 120 121 122 /*********************************************************************/ 123 124 /**************/ 125 /*** ALGO E ***/ 126 /**************/ 127 128 /* similar to shift_unroll_cc() - but, have fixed and limited precalc size 129 * idea: smaller cache usage by table 130 * size must be multiple of CSDR_SHIFT_LIMITED_SIMD (= 4) 131 */ 132 #define PF_SHIFT_LIMITED_UNROLL_SIZE 128 133 #define PF_SHIFT_LIMITED_SIMD_SZ 4 134 135 typedef struct shift_limited_unroll_data_s 136 { 137 float dcos[PF_SHIFT_LIMITED_UNROLL_SIZE]; 138 float dsin[PF_SHIFT_LIMITED_UNROLL_SIZE]; 139 complexf complex_phase; 140 float phase_increment; 141 } shift_limited_unroll_data_t; 142 143 shift_limited_unroll_data_t shift_limited_unroll_init(float rate); 144 /* size must be multiple of PF_SHIFT_LIMITED_SIMD_SZ */ 145 /* starting_phase for next call is kept internal in state */ 146 void shift_limited_unroll_cc(const complexf *input, complexf* output, int size, shift_limited_unroll_data_t* d); 147 void shift_limited_unroll_inp_c(complexf* in_out, int size, shift_limited_unroll_data_t* d); 148 149 150 /*********************************************************************/ 151 152 /**************/ 153 /*** ALGO F ***/ 154 /**************/ 155 156 typedef struct shift_limited_unroll_A_sse_data_s 157 { 158 /* small/limited trig table */ 159 float dcos[PF_SHIFT_LIMITED_UNROLL_SIZE+PF_SHIFT_LIMITED_SIMD_SZ]; 160 float dsin[PF_SHIFT_LIMITED_UNROLL_SIZE+PF_SHIFT_LIMITED_SIMD_SZ]; 161 /* 4 times complex phase */ 162 float phase_state_i[PF_SHIFT_LIMITED_SIMD_SZ]; 163 float phase_state_q[PF_SHIFT_LIMITED_SIMD_SZ]; 164 /* N_cplx_per_block times increment - for future parallel variants */ 165 float dcos_blk; 166 float dsin_blk; 167 /* */ 168 float phase_increment; 169 } shift_limited_unroll_A_sse_data_t; 170 171 shift_limited_unroll_A_sse_data_t shift_limited_unroll_A_sse_init(float relative_freq, float phase_start_rad); 172 void shift_limited_unroll_A_sse_inp_c(complexf* in_out, int N_cplx, shift_limited_unroll_A_sse_data_t* d); 173 174 175 /*********************************************************************/ 176 177 /**************/ 178 /*** ALGO G ***/ 179 /**************/ 180 181 typedef struct shift_limited_unroll_B_sse_data_s 182 { 183 /* small/limited trig table */ 184 float dtrig[PF_SHIFT_LIMITED_UNROLL_SIZE+PF_SHIFT_LIMITED_SIMD_SZ]; 185 /* 4 times complex phase */ 186 float phase_state_i[PF_SHIFT_LIMITED_SIMD_SZ]; 187 float phase_state_q[PF_SHIFT_LIMITED_SIMD_SZ]; 188 /* N_cplx_per_block times increment - for future parallel variants */ 189 float dcos_blk; 190 float dsin_blk; 191 /* */ 192 float phase_increment; 193 } shift_limited_unroll_B_sse_data_t; 194 195 shift_limited_unroll_B_sse_data_t shift_limited_unroll_B_sse_init(float relative_freq, float phase_start_rad); 196 void shift_limited_unroll_B_sse_inp_c(complexf* in_out, int N_cplx, shift_limited_unroll_B_sse_data_t* d); 197 198 /*********************************************************************/ 199 200 /**************/ 201 /*** ALGO H ***/ 202 /**************/ 203 204 typedef struct shift_limited_unroll_C_sse_data_s 205 { 206 /* small/limited trig table - interleaved: 4 cos, 4 sin, 4 cos, .. */ 207 float dinterl_trig[2*(PF_SHIFT_LIMITED_UNROLL_SIZE+PF_SHIFT_LIMITED_SIMD_SZ)]; 208 /* 4 times complex phase */ 209 float phase_state_i[PF_SHIFT_LIMITED_SIMD_SZ]; 210 float phase_state_q[PF_SHIFT_LIMITED_SIMD_SZ]; 211 /* N_cplx_per_block times increment - for future parallel variants */ 212 float dcos_blk; 213 float dsin_blk; 214 /* */ 215 float phase_increment; 216 } shift_limited_unroll_C_sse_data_t; 217 218 shift_limited_unroll_C_sse_data_t shift_limited_unroll_C_sse_init(float relative_freq, float phase_start_rad); 219 void shift_limited_unroll_C_sse_inp_c(complexf* in_out, int N_cplx, shift_limited_unroll_C_sse_data_t* d); 220 221 222 223 /*********************************************************************/ 224 225 /**************/ 226 /*** ALGO I ***/ 227 /**************/ 228 229 /* Recursive Quadrature Oscillator functions "recursive_osc" 230 * see https://www.vicanek.de/articles/QuadOsc.pdf 231 */ 232 #define PF_SHIFT_RECURSIVE_SIMD_SZ 8 233 typedef struct shift_recursive_osc_s 234 { 235 float u_cos[PF_SHIFT_RECURSIVE_SIMD_SZ]; 236 float v_sin[PF_SHIFT_RECURSIVE_SIMD_SZ]; 237 } shift_recursive_osc_t; 238 239 typedef struct shift_recursive_osc_conf_s 240 { 241 float k1; 242 float k2; 243 } shift_recursive_osc_conf_t; 244 245 void shift_recursive_osc_init(float rate, float starting_phase, shift_recursive_osc_conf_t *conf, shift_recursive_osc_t *state); 246 void shift_recursive_osc_update_rate(float rate, shift_recursive_osc_conf_t *conf, shift_recursive_osc_t* state); 247 248 /* size must be multiple of PF_SHIFT_LIMITED_SIMD_SZ */ 249 /* starting_phase for next call is kept internal in state */ 250 void shift_recursive_osc_cc(const complexf *input, complexf* output, int size, const shift_recursive_osc_conf_t *conf, shift_recursive_osc_t* state); 251 void shift_recursive_osc_inp_c(complexf* output, int size, const shift_recursive_osc_conf_t *conf, shift_recursive_osc_t* state); 252 void gen_recursive_osc_c(complexf* output, int size, const shift_recursive_osc_conf_t *conf, shift_recursive_osc_t* state); 253 254 /*********************************************************************/ 255 256 /**************/ 257 /*** ALGO J ***/ 258 /**************/ 259 260 #define PF_SHIFT_RECURSIVE_SIMD_SSE_SZ 4 261 typedef struct shift_recursive_osc_sse_s 262 { 263 float u_cos[PF_SHIFT_RECURSIVE_SIMD_SSE_SZ]; 264 float v_sin[PF_SHIFT_RECURSIVE_SIMD_SSE_SZ]; 265 } shift_recursive_osc_sse_t; 266 267 typedef struct shift_recursive_osc_sse_conf_s 268 { 269 float k1; 270 float k2; 271 } shift_recursive_osc_sse_conf_t; 272 273 void shift_recursive_osc_sse_init(float rate, float starting_phase, shift_recursive_osc_sse_conf_t *conf, shift_recursive_osc_sse_t *state); 274 void shift_recursive_osc_sse_update_rate(float rate, shift_recursive_osc_sse_conf_t *conf, shift_recursive_osc_sse_t* state); 275 void shift_recursive_osc_sse_inp_c(complexf* in_out, int N_cplx, const shift_recursive_osc_sse_conf_t *conf, shift_recursive_osc_sse_t* state_ext); 276 277 278 #ifdef __cplusplus 279 } 280 #endif 281 282