1 /* 2 * Copyright (C) 2020 Linux Studio Plugins Project <https://lsp-plug.in/> 3 * (C) 2020 Vladimir Sadovnikov <sadko4u@gmail.com> 4 * 5 * This file is part of lsp-plugins 6 * Created on: 15 февр. 2017 г. 7 * 8 * lsp-plugins is free software: you can redistribute it and/or modify 9 * it under the terms of the GNU Lesser General Public License as published by 10 * the Free Software Foundation, either version 3 of the License, or 11 * any later version. 12 * 13 * lsp-plugins is distributed in the hope that it will be useful, 14 * but WITHOUT ANY WARRANTY; without even the implied warranty of 15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 16 * GNU Lesser General Public License for more details. 17 * 18 * You should have received a copy of the GNU Lesser General Public License 19 * along with lsp-plugins. If not, see <https://www.gnu.org/licenses/>. 20 */ 21 22 #ifndef DSP_ARCH_NATIVE_COMPLEX_H_ 23 #define DSP_ARCH_NATIVE_COMPLEX_H_ 24 25 #ifndef __DSP_NATIVE_IMPL 26 #error "This header should not be included directly" 27 #endif /* __DSP_NATIVE_IMPL */ 28 29 namespace native 30 { complex_mul2(float * dst_re,float * dst_im,const float * src_re,const float * src_im,size_t count)31 void complex_mul2(float *dst_re, float *dst_im, const float *src_re, const float *src_im, size_t count) 32 { 33 for (size_t i=0; i<count; ++i) 34 { 35 float re = (dst_re[i]) * (src_re[i]) - (dst_im[i]) * (src_im[i]); 36 float im = (dst_re[i]) * (src_im[i]) + (dst_im[i]) * (src_re[i]); 37 dst_re[i] = re; 38 dst_im[i] = im; 39 } 40 } 41 complex_mul3(float * dst_re,float * dst_im,const float * src1_re,const float * src1_im,const float * src2_re,const float * src2_im,size_t count)42 void complex_mul3(float *dst_re, float *dst_im, const float *src1_re, const float *src1_im, const float *src2_re, const float *src2_im, size_t count) 43 { 44 for (size_t i=0; i<count; ++i) 45 { 46 float re = (src1_re[i]) * (src2_re[i]) - (src1_im[i]) * (src2_im[i]); 47 float im = (src1_re[i]) * (src2_im[i]) + (src1_im[i]) * (src2_re[i]); 48 dst_re[i] = re; 49 dst_im[i] = im; 50 } 51 } 52 complex_rcp1(float * dst_re,float * dst_im,size_t count)53 void complex_rcp1(float *dst_re, float *dst_im, size_t count) 54 { 55 while (count--) 56 { 57 float re = *dst_re; 58 float im = *dst_im; 59 float mag = 1.0f / (re * re + im * im); 60 61 *(dst_re++) = re * mag; 62 *(dst_im++) = -im * mag; 63 } 64 } 65 complex_rcp2(float * dst_re,float * dst_im,const float * src_re,const float * src_im,size_t count)66 void complex_rcp2(float *dst_re, float *dst_im, const float *src_re, const float *src_im, size_t count) 67 { 68 while (count--) 69 { 70 float re = *(src_re++); 71 float im = *(src_im++); 72 float mag = 1.0f / (re * re + im * im); 73 74 *(dst_re++) = re * mag; 75 *(dst_im++) = -im * mag; 76 } 77 } 78 complex_cvt2modarg(float * dst_mod,float * dst_arg,const float * src_re,const float * src_im,size_t count)79 void complex_cvt2modarg(float *dst_mod, float *dst_arg, const float *src_re, const float *src_im, size_t count) 80 { 81 while (count--) 82 { 83 float r = *(src_re++); 84 float i = *(src_im++); 85 float r2 = r * r; 86 float i2 = i * i; 87 88 float m = sqrtf(r2 + i2); 89 float a = (i != 0.0f) ? 2.0f * atanf((m - r)/ i) : 90 (r == 0.0f) ? NAN : 91 (r < 0.0f) ? M_PI : 0.0f; 92 93 *(dst_mod++) = m; 94 *(dst_arg++) = a; 95 } 96 } 97 complex_arg(float * dst,const float * re,const float * im,size_t count)98 void complex_arg(float *dst, const float *re, const float *im, size_t count) 99 { 100 while (count--) 101 { 102 float r = *(re++); 103 float i = *(im++); 104 float r2 = r * r; 105 float i2 = i * i; 106 107 float m = sqrtf(r2 + i2); 108 float a = (i != 0.0f) ? 2.0f * atanf((m - r)/ i) : 109 (r == 0.0f) ? NAN : 110 (r < 0.0f) ? M_PI : 0.0f; 111 112 *(dst++) = a; 113 } 114 } 115 complex_cvt2reim(float * dst_re,float * dst_im,const float * src_mod,const float * src_arg,size_t count)116 void complex_cvt2reim(float *dst_re, float *dst_im, const float *src_mod, const float *src_arg, size_t count) 117 { 118 while (count--) 119 { 120 float mod = *(src_mod++); 121 float arg = *(src_arg++); 122 *(dst_re++) = mod * cosf(arg); 123 *(dst_im++) = mod * sinf(arg); 124 } 125 } 126 complex_mod(float * dst_mod,const float * src_re,const float * src_im,size_t count)127 void complex_mod(float *dst_mod, const float *src_re, const float *src_im, size_t count) 128 { 129 while (count--) 130 { 131 float re = *(src_re++); 132 float im = *(src_im++); 133 *(dst_mod++) = sqrtf(re*re + im*im); 134 } 135 } 136 complex_div2(float * dst_re,float * dst_im,const float * src_re,const float * src_im,size_t count)137 void complex_div2(float *dst_re, float *dst_im, const float *src_re, const float *src_im, size_t count) 138 { 139 for (size_t i=0; i<count; ++i) 140 { 141 float re = src_re[i] * dst_re[i] + src_im[i] * dst_im[i]; 142 float im = src_re[i] * dst_im[i] + src_im[i] * dst_re[i]; 143 float n = 1.0f / (src_re[i] * src_re[i] + src_im[i] * src_im[i]); 144 145 dst_re[i] = re * n; 146 dst_im[i] = -im * n; 147 } 148 } 149 complex_rdiv2(float * dst_re,float * dst_im,const float * src_re,const float * src_im,size_t count)150 void complex_rdiv2(float *dst_re, float *dst_im, const float *src_re, const float *src_im, size_t count) 151 { 152 for (size_t i=0; i<count; ++i) 153 { 154 float re = src_re[i] * dst_re[i] + src_im[i] * dst_im[i]; 155 float im = src_re[i] * dst_im[i] + src_im[i] * dst_re[i]; 156 float n = 1.0f / (dst_re[i] * dst_re[i] + dst_im[i] * dst_im[i]); 157 158 dst_re[i] = re * n; 159 dst_im[i] = -im * n; 160 } 161 } 162 complex_div3(float * dst_re,float * dst_im,const float * t_re,const float * t_im,const float * b_re,const float * b_im,size_t count)163 void complex_div3(float *dst_re, float *dst_im, const float *t_re, const float *t_im, const float *b_re, const float *b_im, size_t count) 164 { 165 for (size_t i=0; i<count; ++i) 166 { 167 float re = t_re[i] * b_re[i] + t_im[i] * b_im[i]; 168 float im = t_re[i] * b_im[i] + t_im[i] * b_re[i]; 169 float n = 1.0f / (b_re[i] * b_re[i] + b_im[i] * b_im[i]); 170 171 dst_re[i] = re * n; 172 dst_im[i] = -im * n; 173 } 174 } 175 176 } 177 178 #endif /* DSP_ARCH_NATIVE_COMPLEX_H_ */ 179