GNU Radio 3.4.0 C++ API
volk_32fc_x2_conjugate_dot_prod_32fc_a16.h
Go to the documentation of this file.
00001 #ifndef INCLUDED_volk_32fc_x2_conjugate_dot_prod_32fc_a16_H
00002 #define INCLUDED_volk_32fc_x2_conjugate_dot_prod_32fc_a16_H
00003 
00004 #include<volk/volk_complex.h>
00005 #include<stdio.h>
00006 
00007 
00008 #if LV_HAVE_GENERIC
00009 
00010 
00011 static inline void volk_32fc_x2_conjugate_dot_prod_32fc_a16_generic(lv_32fc_t* result, const lv_32fc_t* input, const lv_32fc_t* taps, unsigned int num_bytes) {
00012   
00013   float * res = (float*) result;
00014   float * in = (float*) input;
00015   float * tp = (float*) taps;
00016   unsigned int n_2_ccomplex_blocks = num_bytes >> 4;
00017   unsigned int isodd = (num_bytes >> 3) &1;
00018   
00019   
00020   
00021   float sum0[2] = {0,0};
00022   float sum1[2] = {0,0};
00023   int i = 0;
00024 
00025   
00026   for(i = 0; i < n_2_ccomplex_blocks; ++i) {
00027     
00028 
00029     sum0[0] += in[0] * tp[0] + in[1] * tp[1];
00030     sum0[1] += (-in[0] * tp[1]) + in[1] * tp[0];
00031     sum1[0] += in[2] * tp[2] + in[3] * tp[3];
00032     sum1[1] += (-in[2] * tp[3]) + in[3] * tp[2];
00033     
00034     
00035     in += 4;
00036     tp += 4;
00037 
00038   }
00039  
00040   
00041   res[0] = sum0[0] + sum1[0];
00042   res[1] = sum0[1] + sum1[1];
00043   
00044   
00045   
00046   for(i = 0; i < isodd; ++i) {
00047 
00048 
00049     *result += input[(num_bytes >> 3) - 1] * lv_conj(taps[(num_bytes >> 3) - 1]);
00050 
00051   }
00052   /*
00053   for(i = 0; i < num_bytes >> 3; ++i) {
00054     *result += input[i] * conjf(taps[i]);
00055   }
00056   */
00057 }
00058 
00059 #endif /*LV_HAVE_GENERIC*/
00060 
00061 
00062 #if LV_HAVE_SSE && LV_HAVE_64
00063 
00064 
00065 static inline void volk_32fc_x2_conjugate_dot_prod_32fc_a16_sse(lv_32fc_t* result, const lv_32fc_t* input, const lv_32fc_t* taps, unsigned int num_bytes) {
00066   
00067   static const uint32_t conjugator[4] __attribute__((aligned(16)))= {0x00000000, 0x80000000, 0x00000000, 0x80000000};
00068   
00069 
00070 
00071 
00072   asm volatile 
00073     (
00074      "#  ccomplex_conjugate_dotprod_generic (float* result, const float *input,\n\t"
00075      "#                         const float *taps, unsigned num_bytes)\n\t"
00076      "#    float sum0 = 0;\n\t"
00077      "#    float sum1 = 0;\n\t"
00078      "#    float sum2 = 0;\n\t"
00079      "#    float sum3 = 0;\n\t"
00080      "#    do {\n\t"
00081      "#      sum0 += input[0] * taps[0] - input[1] * taps[1];\n\t"
00082      "#      sum1 += input[0] * taps[1] + input[1] * taps[0];\n\t"
00083      "#      sum2 += input[2] * taps[2] - input[3] * taps[3];\n\t"
00084      "#      sum3 += input[2] * taps[3] + input[3] * taps[2];\n\t"
00085      "#      input += 4;\n\t"
00086      "#      taps += 4;  \n\t"
00087      "#    } while (--n_2_ccomplex_blocks != 0);\n\t"
00088      "#    result[0] = sum0 + sum2;\n\t"
00089      "#    result[1] = sum1 + sum3;\n\t"
00090      "# TODO: prefetch and better scheduling\n\t"
00091      "  xor    %%r9,  %%r9\n\t"
00092      "  xor    %%r10, %%r10\n\t"
00093      "  movq   %[conjugator], %%r9\n\t"
00094      "  movq   %%rcx, %%rax\n\t"
00095      "  movaps 0(%%r9), %%xmm8\n\t"
00096      "  movq   %%rcx, %%r8\n\t"
00097      "  movq   %[rsi],  %%r9\n\t"
00098      "  movq   %[rdx], %%r10\n\t"
00099      "  xorps   %%xmm6, %%xmm6          # zero accumulators\n\t"
00100      "  movaps  0(%%r9), %%xmm0\n\t"
00101      "  xorps   %%xmm7, %%xmm7          # zero accumulators\n\t"
00102      "  movups  0(%%r10), %%xmm2\n\t"
00103      "  shr     $5, %%rax               # rax = n_2_ccomplex_blocks / 2\n\t"
00104      "  shr     $4, %%r8\n\t"
00105      "  xorps  %%xmm8, %%xmm2\n\t"
00106      "  jmp     .%=L1_test\n\t"
00107      "  # 4 taps / loop\n\t"
00108      "  # something like ?? cycles / loop\n\t"
00109      ".%=Loop1: \n\t"
00110      "# complex prod: C += A * B,  w/ temp Z & Y (or B), xmmPN=$0x8000000080000000\n\t"
00111      "# movaps  (%%r9), %%xmmA\n\t"
00112      "# movaps  (%%r10), %%xmmB\n\t"
00113      "# movaps  %%xmmA, %%xmmZ\n\t"
00114      "# shufps  $0xb1, %%xmmZ, %%xmmZ   # swap internals\n\t"
00115      "# mulps   %%xmmB, %%xmmA\n\t"
00116      "# mulps   %%xmmZ, %%xmmB\n\t"
00117      "# # SSE replacement for: pfpnacc %%xmmB, %%xmmA\n\t"
00118      "# xorps   %%xmmPN, %%xmmA\n\t"
00119      "# movaps  %%xmmA, %%xmmZ\n\t"
00120      "# unpcklps %%xmmB, %%xmmA\n\t"
00121      "# unpckhps %%xmmB, %%xmmZ\n\t"
00122      "# movaps  %%xmmZ, %%xmmY\n\t"
00123      "# shufps  $0x44, %%xmmA, %%xmmZ   # b01000100\n\t"
00124      "# shufps  $0xee, %%xmmY, %%xmmA   # b11101110\n\t"
00125      "# addps   %%xmmZ, %%xmmA\n\t"
00126      "# addps   %%xmmA, %%xmmC\n\t"
00127      "# A=xmm0, B=xmm2, Z=xmm4\n\t"
00128      "# A'=xmm1, B'=xmm3, Z'=xmm5\n\t"
00129      "  movaps  16(%%r9), %%xmm1\n\t"
00130      "  movaps  %%xmm0, %%xmm4\n\t"
00131      "  mulps   %%xmm2, %%xmm0\n\t"
00132      "  shufps  $0xb1, %%xmm4, %%xmm4   # swap internals\n\t"
00133      "  movaps  16(%%r10), %%xmm3\n\t"
00134      "  movaps  %%xmm1, %%xmm5\n\t"
00135      "  xorps   %%xmm8, %%xmm3\n\t"
00136      "  addps   %%xmm0, %%xmm6\n\t"
00137      "  mulps   %%xmm3, %%xmm1\n\t"
00138      "  shufps  $0xb1, %%xmm5, %%xmm5   # swap internals\n\t"
00139      "  addps   %%xmm1, %%xmm6\n\t"
00140      "  mulps   %%xmm4, %%xmm2\n\t"
00141      "  movaps  32(%%r9), %%xmm0\n\t"
00142      "  addps   %%xmm2, %%xmm7\n\t"
00143      "  mulps   %%xmm5, %%xmm3\n\t"
00144      "  add     $32, %%r9\n\t"
00145      "  movaps  32(%%r10), %%xmm2\n\t"
00146      "  addps   %%xmm3, %%xmm7\n\t"
00147      "  add     $32, %%r10\n\t"
00148      "  xorps   %%xmm8, %%xmm2\n\t"
00149      ".%=L1_test:\n\t"
00150      "  dec     %%rax\n\t"
00151      "  jge     .%=Loop1\n\t"
00152      "  # We've handled the bulk of multiplies up to here.\n\t"
00153      "  # Let's sse if original n_2_ccomplex_blocks was odd.\n\t"
00154      "  # If so, we've got 2 more taps to do.\n\t"
00155      "  and     $1, %%r8\n\t"
00156      "  je      .%=Leven\n\t"
00157      "  # The count was odd, do 2 more taps.\n\t"
00158      "  # Note that we've already got mm0/mm2 preloaded\n\t"
00159      "  # from the main loop.\n\t"
00160      "  movaps  %%xmm0, %%xmm4\n\t"
00161      "  mulps   %%xmm2, %%xmm0\n\t"
00162      "  shufps  $0xb1, %%xmm4, %%xmm4   # swap internals\n\t"
00163      "  addps   %%xmm0, %%xmm6\n\t"
00164      "  mulps   %%xmm4, %%xmm2\n\t"
00165      "  addps   %%xmm2, %%xmm7\n\t"
00166      ".%=Leven:\n\t"
00167      "  # neg inversor\n\t"
00168      "  xorps   %%xmm1, %%xmm1\n\t"
00169      "  mov     $0x80000000, %%r9\n\t"
00170      "  movd    %%r9, %%xmm1\n\t"
00171      "  shufps  $0x11, %%xmm1, %%xmm1   # b00010001 # 0 -0 0 -0\n\t"
00172      "  # pfpnacc\n\t"
00173      "  xorps   %%xmm1, %%xmm6\n\t"
00174      "  movaps  %%xmm6, %%xmm2\n\t"
00175      "  unpcklps %%xmm7, %%xmm6\n\t"
00176      "  unpckhps %%xmm7, %%xmm2\n\t"
00177      "  movaps  %%xmm2, %%xmm3\n\t"
00178      "  shufps  $0x44, %%xmm6, %%xmm2   # b01000100\n\t"
00179      "  shufps  $0xee, %%xmm3, %%xmm6   # b11101110\n\t"
00180      "  addps   %%xmm2, %%xmm6\n\t"
00181      "                                  # xmm6 = r1 i2 r3 i4\n\t"
00182      "  movhlps %%xmm6, %%xmm4          # xmm4 = r3 i4 ?? ??\n\t"
00183      "  addps   %%xmm4, %%xmm6          # xmm6 = r1+r3 i2+i4 ?? ??\n\t"
00184      "  movlps  %%xmm6, (%[rdi])                # store low 2x32 bits (complex) to memory\n\t"
00185      :
00186      :[rsi] "r" (input), [rdx] "r" (taps), "c" (num_bytes), [rdi] "r" (result), [conjugator] "r" (conjugator)
00187      :"rax", "r8", "r9", "r10"
00188      );
00189   
00190   
00191   int getem = num_bytes % 16;
00192   
00193   
00194   for(; getem > 0; getem -= 8) {
00195   
00196     
00197     *result += (input[(num_bytes >> 3) - 1] * lv_conj(taps[(num_bytes >> 3) - 1]));
00198   
00199   }
00200 
00201   return;
00202 }  
00203 #endif
00204 
00205 #if LV_HAVE_SSE && LV_HAVE_32
00206 static inline void volk_32fc_x2_conjugate_dot_prod_32fc_a16_sse_32(lv_32fc_t* result, const lv_32fc_t* input, const lv_32fc_t* taps, unsigned int num_bytes) {
00207   
00208   static const uint32_t conjugator[4] __attribute__((aligned(16)))= {0x00000000, 0x80000000, 0x00000000, 0x80000000};
00209 
00210   int bound = num_bytes >> 4;
00211   int leftovers = num_bytes % 16;
00212 
00213   
00214   asm volatile 
00215     (
00216      "  #pushl  %%ebp\n\t"
00217      "  #movl   %%esp, %%ebp\n\t"
00218      "  #movl   12(%%ebp), %%eax                # input\n\t"
00219      "  #movl   16(%%ebp), %%edx                # taps\n\t"
00220      "  #movl   20(%%ebp), %%ecx                # n_bytes\n\t"
00221      "  movaps  0(%[conjugator]), %%xmm1\n\t"
00222      "  xorps   %%xmm6, %%xmm6          # zero accumulators\n\t"
00223      "  movaps  0(%[eax]), %%xmm0\n\t"
00224      "  xorps   %%xmm7, %%xmm7          # zero accumulators\n\t"
00225      "  movaps  0(%[edx]), %%xmm2\n\t"
00226      "  movl    %[ecx], (%[out])\n\t"
00227      "  shrl    $5, %[ecx]              # ecx = n_2_ccomplex_blocks / 2\n\t"
00228      
00229      "  xorps   %%xmm1, %%xmm2\n\t"
00230      "  jmp     .%=L1_test\n\t"
00231      "  # 4 taps / loop\n\t"
00232      "  # something like ?? cycles / loop\n\t"
00233      ".%=Loop1: \n\t"
00234      "# complex prod: C += A * B,  w/ temp Z & Y (or B), xmmPN=$0x8000000080000000\n\t"
00235      "# movaps  (%[eax]), %%xmmA\n\t"
00236      "# movaps  (%[edx]), %%xmmB\n\t"
00237      "# movaps  %%xmmA, %%xmmZ\n\t"
00238      "# shufps  $0xb1, %%xmmZ, %%xmmZ   # swap internals\n\t"
00239      "# mulps   %%xmmB, %%xmmA\n\t"
00240      "# mulps   %%xmmZ, %%xmmB\n\t"
00241      "# # SSE replacement for: pfpnacc %%xmmB, %%xmmA\n\t"
00242      "# xorps   %%xmmPN, %%xmmA\n\t"
00243      "# movaps  %%xmmA, %%xmmZ\n\t"
00244      "# unpcklps %%xmmB, %%xmmA\n\t"
00245      "# unpckhps %%xmmB, %%xmmZ\n\t"
00246      "# movaps  %%xmmZ, %%xmmY\n\t"
00247      "# shufps  $0x44, %%xmmA, %%xmmZ   # b01000100\n\t"
00248      "# shufps  $0xee, %%xmmY, %%xmmA   # b11101110\n\t"
00249      "# addps   %%xmmZ, %%xmmA\n\t"
00250      "# addps   %%xmmA, %%xmmC\n\t"
00251      "# A=xmm0, B=xmm2, Z=xmm4\n\t"
00252      "# A'=xmm1, B'=xmm3, Z'=xmm5\n\t"
00253      "  movaps  16(%[edx]), %%xmm3\n\t"
00254      "  movaps  %%xmm0, %%xmm4\n\t"
00255      "  xorps   %%xmm1, %%xmm3\n\t"
00256      "  mulps   %%xmm2, %%xmm0\n\t"
00257      "  movaps  16(%[eax]), %%xmm1\n\t"
00258      "  shufps  $0xb1, %%xmm4, %%xmm4   # swap internals\n\t"
00259      "  movaps  %%xmm1, %%xmm5\n\t"
00260      "  addps   %%xmm0, %%xmm6\n\t"
00261      "  mulps   %%xmm3, %%xmm1\n\t"
00262      "  shufps  $0xb1, %%xmm5, %%xmm5   # swap internals\n\t"
00263      "  addps   %%xmm1, %%xmm6\n\t"
00264      "  movaps  0(%[conjugator]), %%xmm1\n\t"
00265      "  mulps   %%xmm4, %%xmm2\n\t"
00266      "  movaps  32(%[eax]), %%xmm0\n\t"
00267      "  addps   %%xmm2, %%xmm7\n\t"
00268      "  mulps   %%xmm5, %%xmm3\n\t"
00269      "  addl    $32, %[eax]\n\t"
00270      "  movaps  32(%[edx]), %%xmm2\n\t"
00271      "  addps   %%xmm3, %%xmm7\n\t"
00272      "  xorps   %%xmm1, %%xmm2\n\t"
00273      "  addl    $32, %[edx]\n\t"
00274      ".%=L1_test:\n\t"
00275      "  decl    %[ecx]\n\t"
00276      "  jge     .%=Loop1\n\t"
00277      "  # We've handled the bulk of multiplies up to here.\n\t"
00278      "  # Let's sse if original n_2_ccomplex_blocks was odd.\n\t"
00279      "  # If so, we've got 2 more taps to do.\n\t"
00280      "  movl    0(%[out]), %[ecx]               # n_2_ccomplex_blocks\n\t"
00281      "  shrl    $4, %[ecx]\n\t"
00282      "  andl    $1, %[ecx]\n\t"
00283      "  je      .%=Leven\n\t"
00284      "  # The count was odd, do 2 more taps.\n\t"
00285      "  # Note that we've already got mm0/mm2 preloaded\n\t"
00286      "  # from the main loop.\n\t"
00287      "  movaps  %%xmm0, %%xmm4\n\t"
00288      "  mulps   %%xmm2, %%xmm0\n\t"
00289      "  shufps  $0xb1, %%xmm4, %%xmm4   # swap internals\n\t"
00290      "  addps   %%xmm0, %%xmm6\n\t"
00291      "  mulps   %%xmm4, %%xmm2\n\t"
00292      "  addps   %%xmm2, %%xmm7\n\t"
00293      ".%=Leven:\n\t"
00294      "  # neg inversor\n\t"
00295      "  #movl 8(%%ebp), %[eax] \n\t"
00296      "  xorps   %%xmm1, %%xmm1\n\t"
00297      "  movl    $0x80000000, (%[out])\n\t"
00298      "  movss   (%[out]), %%xmm1\n\t"
00299      "  shufps  $0x11, %%xmm1, %%xmm1   # b00010001 # 0 -0 0 -0\n\t"
00300      "  # pfpnacc\n\t"
00301      "  xorps   %%xmm1, %%xmm6\n\t"
00302      "  movaps  %%xmm6, %%xmm2\n\t"
00303      "  unpcklps %%xmm7, %%xmm6\n\t"
00304      "  unpckhps %%xmm7, %%xmm2\n\t"
00305      "  movaps  %%xmm2, %%xmm3\n\t"
00306      "  shufps  $0x44, %%xmm6, %%xmm2   # b01000100\n\t"
00307      "  shufps  $0xee, %%xmm3, %%xmm6   # b11101110\n\t"
00308      "  addps   %%xmm2, %%xmm6\n\t"
00309      "                                  # xmm6 = r1 i2 r3 i4\n\t"
00310      "  #movl   8(%%ebp), %[eax]                # @result\n\t"
00311      "  movhlps %%xmm6, %%xmm4          # xmm4 = r3 i4 ?? ??\n\t"
00312      "  addps   %%xmm4, %%xmm6          # xmm6 = r1+r3 i2+i4 ?? ??\n\t"
00313      "  movlps  %%xmm6, (%[out])                # store low 2x32 bits (complex) to memory\n\t"
00314      "  #popl   %%ebp\n\t"
00315      :
00316      : [eax] "r" (input), [edx] "r" (taps), [ecx] "r" (num_bytes), [out] "r" (result), [conjugator] "r" (conjugator)
00317      );
00318 
00319   
00320   
00321   
00322   printf("%d, %d\n", leftovers, bound);
00323   
00324   for(; leftovers > 0; leftovers -= 8) {
00325     
00326     
00327     *result += (input[(bound << 1)] * lv_conj(taps[(bound << 1)]));
00328     
00329   }
00330   
00331   return;
00332   
00333   
00334   
00335 
00336   
00337   
00338 }
00339 
00340 #endif /*LV_HAVE_SSE*/  
00341 
00342 
00343 
00344 #endif /*INCLUDED_volk_32fc_x2_conjugate_dot_prod_32fc_a16_H*/