开通博客赚积分
发布资源赚积分
分类

源码开发语言/平台
当前位置：首页 > 源码/资料 > 通讯/手机编程 > 语音压缩 > 查看源码
lpc_lib.c
资源名称：melp.rar　[点击查看]
上传用户：csczyc
上传日期：2021-02-19
资源大小：1051k
文件大小：31k
源码类别：

语音压缩

开发平台：

C/C++

lpc_lib.c：源码内容
  2. #include        <stdio.h>
  3. #include        <stdlib.h>
  4. #include        <math.h>
  5. #include "spbstd.h"
  6. #include "mathhalf.h"
  7. #include "mathdp31.h"
  8. #include "math_lib.h"
  9. #include "wmops.h"
  10. #include "mat.h"
  11. #include "lpc.h"
  12. #include "constant.h"
  14. extern Shortword frames;
  15. extern FILE *fp_outdat;
  17. #define PRINT 1
  18. #define ALMOST_ONE_Q14   16382  
  20. Shortword lagw_cof[10] = {
  21. 32756,32721,32663,32582,32478,32351,32201,32030,31837,31622};
  23. void lpc_acor(Shortword input[],Shortword win_cof[],
  24.       Shortword r[],Shortword hf_correction,
  25.       Shortword lag_win,Shortword order,Shortword npts)
  26. {
  27.     Shortword i,j;
  28.     Longword L_temp;
  29.     Shortword *inputw;
  30.     Shortword norm_var, scale_fact, temp;
  32. Longword L_temp1,L_temp2,L_Prod;
  33. Shortword temp1;
  35.     MEM_ALLOC(MALLOC,inputw,npts,Shortword);
  37.     if (win_cof) {
  38.       for(i=0; i < npts; i++) {
  39.    L_temp1 = _shr2((Longword)(input[i]),(Longword)4);
  40.    L_temp1 = _smpy(win_cof[i],(Shortword)(0x0000ffffL & L_temp1));
  41.    inputw[i] = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  42.       }
  43.     }
  44.     else {
  45.       v_equ_shr(inputw,input,4,npts);
  46.     }    
  47.    
  48.     L_temp = L_v_magsq(inputw,npts,0,1);
  49.     if (L_temp) {   
  50.       L_temp1 = _norm(L_temp);
  51.   if (L_temp == 0){
  52.   L_temp1 = (Shortword)0;
  53.   } 
  54.   L_temp1 = _shr2(L_temp1,1);
  55.   norm_var = (Shortword)(0x0000ffffL & L_temp1);
  56.   norm_var = sub(4,norm_var);
  58.       if (norm_var < 0)
  59. norm_var = 0;
  60.     }
  61.     else
  62.       norm_var = 0;   
  63.     L_temp = L_v_magsq(inputw,npts,0,1);
  64.     if (L_temp > 0) {
  65.     
  66.       L_temp1 = _norm(L_temp);
  67.   norm_var = (Shortword) (0x0000ffffL & L_temp1);
  68.       norm_var = sub(norm_var,1);
  69.       L_temp = _sshvl(L_temp,norm_var);      
  70.       L_temp = _sadd(L_temp,L_mpy_ls(L_temp,hf_correction));
  71.      
  72.       temp1 = (Shortword) (0x0000ffffL & (L_temp >> 16));
  73.   L_temp1 = (Longword) temp1 << 16;
  74.   L_temp2 = _norm(L_temp1);
  75.   temp = (Shortword) (0x0000ffffL & L_temp2);
  76.   if (L_temp1 == 0){
  77.   temp = (Shortword)0;
  78.   }
  79.   L_temp = _sshvl(L_temp,temp);
  80.   L_temp1 = _sadd2((Longword)norm_var,(Longword)temp);
  81.       norm_var = (Shortword) (0x0000ffffL & L_temp1);
  82.   L_Prod = _sadd(L_temp, 0x00008000L); 
  83.       r[0] = (Shortword)(0x0000ffffL & (L_Prod >> 16));      
  84.       scale_fact = divide_s(ALMOST_ONE_Q14,r[0]);
  85.       L_temp = _sshvl(L_mpy_ls(L_temp,scale_fact),1);
  86.   L_Prod = _sadd(L_temp, 0x00008000L); /* round MSP */
  87.       r[0] = (Shortword)(0x0000ffffL & (L_Prod >> 16));
  89.     }
  90.     else {
  91.       norm_var = 0;   
  93.       r[0] = ONE_Q15;  
  94.       scale_fact = 0;      
  96.     }    
  97.     for(j=1; j <= order; j++)
  98.     {
  99. L_temp = 0; 
  100.         for(i=j; i < npts; i++) {
  101.     L_temp = _sadd(L_temp,_smpy(inputw[i],inputw[i-j]));
  102. }
  103. L_temp = _sshvl(L_temp,norm_var);
  105. L_temp = _sshvl(L_mpy_ls(L_temp,scale_fact),1);
  107. L_temp = L_mpy_ls(L_temp,lagw_cof[j-1]);
  108. L_Prod = _sadd(L_temp, 0x00008000L); 
  109.     r[j] = (Shortword)(0x0000ffffL & (L_Prod >> 16));
  110.     }
  111.     
  112.     MEM_FREE(FREE,inputw);
  114. }
  116. #define LOW_LIMIT 54
  118. Longword lpc_aejw(Shortword *a,Shortword w, Shortword p)
  119. {
  120.     Shortword i;
  121.     Shortword c_re,c_im;
  122.     Shortword cs,sn,temp;
  123.     Shortword temp1,temp2;
  124.     Longword L_temp;
  125. Longword L_temp1;   
  127.     if (p==0)
  128.         return((Longword)ONE_Q19);  
  129.     cs = cos_fxp(w); 
  130.     sn = negate(sin_fxp(w));
  131.     L_temp1 = _smpy(cs,a[p]);
  132.     L_temp1 = _shr2(L_temp1,3);
  133. c_re = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  134. L_temp1 = _smpy(sn,a[p]);
  135.     L_temp1 = _shr2(L_temp1,3);
  136. c_im = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  137.     for(i=p-1; i > 0; i--)
  138.     {
  139.            L_temp1 = _shr2((Longword)a[i],3);
  140.    L_temp1 = _sadd2((Longword)c_re,L_temp1);
  141.    c_re = (Shortword) (0x0000ffffL & L_temp1);
  142. temp=c_im;
  143.     L_temp1 = _smpy(cs,temp);
  144. temp1 = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  145. L_temp1 = _smpy(sn,c_re);
  146. temp2 = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  147. L_temp1 = _sadd2((Longword)temp1,(Longword)temp2);
  148. c_im = (Shortword) (0x0000ffffL & L_temp1);
  149.     L_temp1 = _smpy(cs,c_re);
  150. temp1 = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  151. L_temp1 = _smpy(sn,temp);
  152. temp2 = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  153.         c_re = sub(temp1,temp2);
  154.     }
  155.    
  156. L_temp1 = _sadd2((Longword)c_re,(Longword)ONE_Q9);
  157. c_re = (Shortword) (0x0000ffffL & L_temp1);
  158.      L_temp = _sadd(_smpy(c_re,c_re),_smpy(c_im,c_im));
  159.     
  160.     if (L_temp<LOW_LIMIT)
  161.       L_temp=(Longword)LOW_LIMIT;
  162.     
  163.     return(L_temp);
  164. } 
  166. #undef LOW_LIMIT
  170. Shortword lpc_bwex(Shortword *a, Shortword *aw, Shortword gamma, Shortword p)
  171. {
  172.     Shortword i;
  173.     Shortword gk;
  175. Longword L_temp1,L_temp2;
  177.     gk = gamma;
  178.     for(i=1; i <= p; i++) {  
  179.       L_temp1 = _smpy(a[i],gk);
  180.   aw[i] = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  181.   L_temp2 = _smpy(gk,gamma);
  182.   gk = (Shortword) (0x0000ffffL & (L_temp2 >> 16));
  183.     }
  184.     return((Shortword)0);
  185. }
  189. #define MAX_LOOPS 10
  191. Shortword lpc_clmp(Shortword *w, Shortword delta, Shortword p)
  192. {
  193.     Shortword i,j,unsorted;
  194.     Shortword temp,d,step1,step2;
  195. Longword L_temp;
  196.     
  197.     for (j=0,unsorted=TRUE; unsorted && (j < MAX_LOOPS); j++)
  198.     {
  199.         for(i=1,unsorted=FALSE; i < p; i++)
  200.             if (w[i] > w[i+1])
  201.             {
  202.                 temp = w[i+1];
  203.                 w[i+1] = w[i];
  204.                 w[i] = temp;
  205.                 unsorted = TRUE;
  206.             }
  207.     }
  209.     
  210.     if (!unsorted) 
  211.     {
  212.         for(j=0; j < MAX_LOOPS; j++)
  213.         {
  214.             for(i=1; i < p; i++)
  215.             {
  216.         /* increment complexity for if statement */
  217. //         compare_nonzero();mark del
  218. d = sub(w[i+1],w[i]);
  219. //      data_move();
  220.                 if (d < delta)
  221.                 {
  222. //     data_move();mark del
  223. //         data_move();mark del
  224.                     step1 = step2 = shr_a(sub(delta,d),1);
  226.     /* increment complexity for if statement */
  227. //     compare_nonzero();    compare_nonzero(); mark del
  228.                     if (i==1 && (w[i] < delta)) {
  229.       step1 = shr_a(w[i],1);
  230. //           data_move();  mark del
  231.     }
  232.                     else {
  233.       /* increment complexity for if statement */
  234. //       compare_nonzero();mark del
  235.       if (i > 1)
  236.       {
  237.         /* increment complexity for if statement */
  238. //         compare_nonzero();mark del
  239. temp = sub(w[i],w[i-1]);
  240.                         if (temp < delta) {
  241.   step1 = 0;   // data_move();mark del
  242. }
  243.                         else {
  244.   /* increment complexity for if statement */
  245. //   compare_nonzero();    compare_nonzero(); mark del
  246.   if (temp < shl_a(delta,1)) {
  247.     step1 = shr_a(sub(temp,delta),1);
  248.   }
  249. }
  250.       }
  251.                     }
  253.     /* increment complexity for if statement */
  254. //     compare_nonzero();    compare_nonzero();  mark del
  255.                     if (i==(p-1) && (w[i+1] > sub(ONE_Q15,delta))) {
  256.                         step2 = shr_a(sub(ONE_Q15,w[i+1]),1);
  257. //                            data_move();mark del
  258.     }
  259.                     else {
  260.       /* increment complexity for if statement */
  261. //       compare_nonzero();mark del
  262.       if (i < (p-1))
  263.       {
  264.         /* increment complexity for if statement */
  265. //         compare_nonzero();mark del
  266. temp = sub(w[i+2],w[i+1]);
  267.                         if (temp < delta) {
  268.   step2 = 0;
  269. //       data_move();
  270. }
  271.                         else {
  272.   /* increment complexity for if statement */
  273. //   compare_nonzero();mark del
  274.   if (temp < shl_a(delta,1)) {
  275.                             step2 = shr_a(sub(temp,delta),1);
  276.   }
  277. }
  278.       }
  279.                     }
  280.                     w[i] = sub(w[i],step1);
  281. //                        data_move();mark del
  282.                    /*  w[i+1] = add(w[i+1],step2); */
  283.     L_temp = _sadd2((Longword)w[i+1],(Longword)step2);
  284.                     w[i+1] = (Shortword) (0x0000ffffL & L_temp);
  285. //                        data_move();mark del
  286.                     
  287.                 }
  288.             }
  289.         }
  290.     }  
  291.     
  292.     /* Debug: check if the minimum separation rule was met */
  293.     /* temp = mult(32440,delta);    temp = 0.99*delta */
  294.   L_temp = _smpy(32440,delta);
  295.   temp = (Shortword) (0x0000ffffL & (L_temp >> 16));
  296.      for(j=1; j < p; j++)                                                       
  297.        if ((w[j+1]-w[j]) < temp)                                                
  298.            (void)fprintf(stderr,"%s: LSPs not separated by enough (line %d)n", 
  299.                __FILE__,__LINE__);                                              
  300.                                                                                 
  301.      if (unsorted)                                                            
  302.        (void)fprintf(stderr,"%s: Fxp LSPs still unsorted (line %d)n",       
  303.                      __FILE__,__LINE__);                                       
  305.     return((Shortword)0);
  306. } /* LPC_CLMP */
  309. /*  
  310.     Name: lpc_schr- Schur recursion (autocorrelations to refl coef)
  311.     Aliases: lpc_schur
  312.     Description:
  313.         Compute reflection coefficients from autocorrelations
  314.         based on schur recursion.  Will also compute predictor
  315.         parameters by calling lpc_refl2pred(3l) if necessary.
  316.     Inputs:
  317.         r- autocorrelation vector (r[0..p]).
  318.         p- order of lpc filter.
  319.     Outputs:
  320.         a-     predictor parameters    (can be NULL)
  321.         k_tmp- reflection coefficients (can be NULL)
  322.     Returns:
  323.         alphap- the minimum residual energy
  324.     Includes:
  325.         spbstd.h
  326.         lpc.h
  327.     See_Also:
  328.         lpc_refl2pred(3l) in lpc.h or lpc(3l)
  330.     Q values:
  331.     r - Q0
  332.     a - Q12
  333.     k_tmp - Q15
  334. */
  336. Shortword lpc_schr(Shortword *r, Shortword *a, Shortword *k_tmp, Shortword p)
  337. {
  338.     Shortword i,j;
  339.     Shortword shift,alphap,*k;
  340.     Longword L_temp, *y1, *y2;
  341.     Shortword n_temp,n_y2,nr1,nr0;
  343. Longword L_temp1;
  344. Shortword temp;
  345.     
  346.  
  347.     MEM_ALLOC(MALLOC,y1,p+2,Longword);
  348.     MEM_ALLOC(MALLOC,y2,p+2,Longword);
  350.     if (k_tmp == NULL)
  351.   MEM_ALLOC(MALLOC,k,p+1,Shortword);
  352.     else
  353. k = k_tmp;
  355.  /*    nr1 = abs_s(r[1]);         */
  356. /*        data_move();mark del */
  357. /*     nr0 = abs_s(r[0]);         */
  358. //        data_move();mark del
  359.       L_temp1 = _abs2((Longword)r[1]);
  360.       nr1 = (Shortword) (0x0000ffffL & L_temp1);
  361.       L_temp1 = _abs2((Longword)r[0]);    
  362.   nr0 = (Shortword) (0x0000ffffL & L_temp1);
  364.     k[1] = divide_s(nr1,nr0);
  365. //        data_move();mark del
  367.     /* increment complexity for if statement */
  368. //    logic(); compare_zero();mark del
  369.     /* if (((r[1] < 0) && (r[0] < 0)) || ((r[1] > 0) && (r[0] > 0))) */
  370.     if ((r[1] ^ r[0]) >= 0)
  371.     {
  372.       k[1] = negate(k[1]);
  373.  //         data_move();mark del
  374.     }
  375.    /*  alphap = mult(r[0],sub(ONE_Q15,mult(k[1],k[1])));  */
  376.      L_temp1 = _smpy(k[1],k[1]);
  377.  temp = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  378.  temp = sub(ONE_Q15,temp);
  379.  L_temp1 = _smpy(r[0],temp);
  380.  alphap = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  383. //       data_move();mark del
  385.     /* y2[1] = L_deposit_h(r[1]); */
  386.    /*  y2[2] = L_add(L_deposit_h(r[0]),L_mult(k[1],r[1])); */
  388.     y2[1] = (Longword) r[1] << 16 ;
  389. L_temp1 = (Longword) r[0] << 16 ;
  390. y2[2] = _sadd(L_temp1,_smpy(k[1],r[1]));
  393.     for(i=2; i <= p; i++)
  394.     {
  395. /* y1[1] = L_deposit_h(r[i]);  */
  396. /*  L_temp = L_deposit_h(r[i]);    */
  397.       y1[1] = (Longword) r[i] << 16;
  398.   L_temp = (Longword) r[i] << 16;
  400.         for(j=1; j < i; j++)
  401.         {  
  402.     y1[j+1] = _sadd(y2[j], L_mpy_ls(L_temp, k[j])); 
  403. //        L_data_move();mark del
  404.     L_temp = _sadd(L_temp, L_mpy_ls(y2[j], k[j]));
  405.     y2[j] = y1[j];   
  406. //      L_data_move();mark del
  407. }
  409.         /* k[i] = -temp/y2[i]; */
  410.         if (L_temp > y2[i]) {
  411.   for (j = i; j <= p; j++)
  412.     k[j] = 0;
  413. #if (PRINT)
  414.   printf("nERROR: numerator bigger than denominator in lpc_schrn");
  415.   printf("Reflection coefficients used in frame %d:n",frames);
  416.   for (i = 1; i <= p; i++) {
  417.     printf(" %5.8f ",(float)k[i]/(1<<15));
  418.   }
  419.   printf("n");
  420. #endif
  421.   break;
  422. }
  424. /*  shift = norm_l(y2[i]);                             */
  425. /*  n_temp = abs_s(extract_h(L_shl(L_temp,shift)));    */
  426.     L_temp1 = _norm(y2[i]);
  427. if (y2[i] == 0){
  428.   L_temp1 = (Longword)0;
  429. }
  430. shift = (Shortword) (0x0000ffffL & L_temp1);
  431.     L_temp1 = _sshvl(L_temp,shift);
  432. L_temp1 = _abs2(L_temp1);
  433. n_temp = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  435. //     data_move();mark del
  436. /* n_y2 = abs_s(extract_h(L_shl(y2[i],shift))); */
  437. L_temp1 = _sshvl(y2[i],shift);
  438. L_temp1 = _abs2(L_temp1);
  439. n_y2 = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  440. //     data_move();mark del
  442. k[i] = divide_s(n_temp,n_y2);
  443. //     data_move();mark del
  444. /* increment complexity for if statement */
  445. // L_logic(); compare_zero();mark del
  446. if ((L_temp ^ y2[i]) >= 0)
  447. {
  448.     k[i] = negate(k[i]); 
  449. //        data_move();mark del
  450. }
  452. y2[i+1] = _sadd(y2[i], L_mpy_ls(L_temp, k[i]));
  453. //     L_data_move();mark del
  454. y2[i] = y1[i];   
  455. //  L_data_move();mark del
  457. /*  alphap = mult(alphap,sub(ONE_Q15,mult(k[i],k[i])));    */
  458.     L_temp1 = _smpy(k[i],k[i]);
  459. temp = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  460. temp = sub(ONE_Q15,temp);
  461. L_temp1 = _smpy(alphap,temp);
  462. alphap = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  463. //   data_move();mark del
  464.     }
  465.     
  466.     if (a != NULL)
  467.     {
  468. lpc_refl2pred(k,a,p);
  469.     }
  471.     if (k_tmp == NULL) 
  472. MEM_FREE(FREE,k);
  474.     MEM_FREE(FREE,y2);
  475.     MEM_FREE(FREE,y1);
  476.     
  477.     return(alphap);   /* alphap in Q15 */
  478. } /* LPC_SCHR */
  481. /* LPC_REFL2PRED
  482.       get predictor coefficients from the reflection coeffs
  483.    Synopsis: lpc_refl2pred(k,a,p)
  485.       Input:
  486.          k- the reflection coeffs
  487.          p- the predictor order
  488.       Output:
  489.          a- the predictor coefficients
  490.    Reference:  Markel and Gray, Linear Prediction of Speech
  492.    Q values:
  493.    k - Q15
  494.    a - Q12
  495. */
  497. Shortword lpc_refl2pred(Shortword *k,Shortword *a,Shortword p)
  498. {
  499.     Shortword i,j;
  500.     Shortword *a1;
  502. Longword L_temp,swRnd;
  503. Shortword temp;
  505.     MEM_ALLOC(MALLOC,a1,p+1,Shortword);
  507.     for(i=1; i <= p; i++)
  508.     {
  509.         /* refl to a recursion */
  510.        /*   a[i] = shift_r(k[i],-3);       a in Q12   */ 
  511.     L_temp = _shr2((Longword)k[i],2);
  512.     swRnd = L_temp & (Longword)0x1;
  513.     L_temp = _shr2((Longword)k[i],3);          
  514.   L_temp = _sadd2(L_temp,swRnd);       
  515.   a[i] = (Shortword) (0x0000ffffL & L_temp);   
  516. //          data_move();mark del
  517.         for(j=1; j < i; j++) {
  518.   a1[j] = a[j];    
  519. //   data_move();mark del
  520. }
  521.         for(j=1; j < i; j++) {
  522.            /*  a[j] = add(a1[j],mult(k[i],a1[i-j]));  */
  523.    L_temp = _smpy(k[i],a1[i-j]);
  524.    temp = (Shortword) (0x0000ffffL & (L_temp >> 16));
  525.    L_temp = _sadd2((Longword)a1[j],(Longword)temp);
  526.    a[j] = (Shortword) (0x0000ffff & L_temp);
  527. //               data_move();mark del
  528. }
  529.     }
  531.     MEM_FREE(FREE,a1);
  532.     return((Shortword)0);
  533. } /* LPC_REFL2PRED */
  536. /* minimum LSP separation-- a global variable */
  537. Shortword lsp_delta = 0;
  539. /* LPC_PRED2LSP
  540.       get LSP coeffs from the predictor coeffs
  541.       Input:
  542.          a- the predictor coefficients
  543.          p- the predictor order
  544.       Output:
  545.          w- the lsp coefficients
  547.       This function uses a DFT to evaluate the P and Q polynomials,
  548.       and is hard-coded to work only for 10th order LPC.
  550.    Q values:
  551.    a - Q12
  552.    w - Q15
  553. */
  554. static void lsp_to_freq(Shortword lsp[restrict], Shortword w[restrict], Shortword order);
  556. #define LPC_ORD 10
  558. Shortword lpc_pred2lsp(Shortword *restrict a,Shortword *restrict w,Shortword p)
  559. {
  560.     Shortword i,j;
  561.     Shortword p_cof[LPC_ORD+1], q_cof[LPC_ORD+1],
  562.               p_freq[LPC_ORD+1], q_freq[LPC_ORD+1];
  563.     Longword L_p_cof[LPC_ORD+1], L_q_cof[LPC_ORD+1];
  564.     Longword L_ai,L_api,L_temp;
  566. Longword L_temp1,L_Prod;
  568.     /* Generate P' and Q' polynomials   */
  569.     L_p_cof[0] = (Longword)ONE_Q26;
  570.     L_q_cof[0] = (Longword)ONE_Q26;
  571.     for ( i = 1; i <= p; i++ ) {
  572.         /* temp = sub(a[i],a[p+1-i]); */     /* temp in Q12 */
  573.        /*  L_ai = L_shr(L_deposit_h(a[i]),2); */
  574.            L_ai =(Longword)(a[i] << 16);
  575.            L_ai = _sshvr(L_ai,2); 
  576. //             L_data_move();mark del
  577. /* L_api = L_shr(L_deposit_h(a[p+1-i]),2); */
  578.  L_api = (Longword) (a[p+1-i] << 16);
  579.  L_api = _sshvr(L_api,2);
  580. //      L_data_move();  mark del
  581. /*  L_temp = L_sub(L_ai,L_api);        L_temp in Q26 */                     
  582. /*         L_p_cof[i] = L_add(L_temp,L_p_cof[i-1]);       L_p_cof in Q26 */ 
  583.     L_temp = _ssub(L_ai,L_api);
  584. L_p_cof[i] = _sadd(L_temp,L_p_cof[i-1]);
  587. /* temp = add(a[i],a[p+1-i]); */
  588. /*  L_temp = L_add(L_ai,L_api);                                               */
  589. /*         L_q_cof[i] = L_sub(L_temp,L_q_cof[i-1]);       L_q_cof in Q26 */ 
  590.      L_temp = _sadd(L_ai,L_api);
  591.         L_q_cof[i] = _ssub(L_temp,L_q_cof[i-1]);      /* L_q_cof in Q26 */
  592.     }
  594.     /* Convert p_cof and q_cof to short */
  595.     for ( i = 0; i <= p; i++ ) {
  596.      /*  p_cof[i] = round(L_p_cof[i]);      p_cof in Q10 */ 
  597.   L_Prod = _sadd(L_p_cof[i], 0x00008000L); /* round MSP */
  598.       p_cof[i] = (Shortword) (0x0000ffffL & (L_Prod >> 16));
  599.       /* q_cof[i] = round(L_q_cof[i]);      q_cof in Q10 */ 
  600.   L_temp1 = _sadd(L_q_cof[i], 0x00008000L); /* round MSP */
  601.       q_cof[i] = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  602.     }
  604.     /* Find root frequencies of LSP polynomials */
  605.     lsp_to_freq(p_cof,p_freq,p);
  606.     lsp_to_freq(q_cof,q_freq,p);
  608.     /* Combine frequencies into single array    */
  609.     w[0] = 0;
  610.     j = 1;
  611.     for (i = 0; i < (shr_a(p,1)); i++) {
  612.         w[j++] = q_freq[i];      
  613. //        data_move();mark del
  614.         w[j++] = p_freq[i];   
  615. //           data_move();mark del
  616.     }
  618.     return((Shortword)0);
  620. } /* LPC_PRED2LSP */
  622. /*                                                              */
  623. /*      Subroutine LSP_TO_FREQ: Calculate Line Spectrum Pair    */
  624. /*      root frequencies from LSP polynomial.                   */
  625. /*                                                              */
  626. /* Q values:
  627.    lsp - Q10
  628.    w - Q15 
  629. */
  630. #define DFTLENGTH       512
  631. #define DFTLENGTH_D2    (DFTLENGTH/2)
  632. #define DFTLENGTH_D4    (DFTLENGTH/4)
  634. static void lsp_to_freq(Shortword lsp[restrict], Shortword w[restrict], Shortword p)
  636. {
  637.     Shortword i;
  638.     Shortword p2, count, prev_less, j;
  639.     static Shortword firstcall = 1;
  640.     Longword mag0, mag1, mag2;
  641.     Shortword s_mag0, s_mag1, s_mag2;
  642.     Shortword *p_cos, *p_endlsp, *p_lsp, *p_begcos, *p_endcos;
  643.     Shortword temp1, temp2;
  644.     Longword L_temp1, L_temp2;
  646. Longword L_temp;
  647.     static Shortword lsp_cos[DFTLENGTH];    /* cosine table */
  648.     static Shortword default_w, default_w0;
  649.     static Shortword cos_tab[129] =
  650.     { 
  651.         32767, 32766, 32758, 32746, 32729, 32706, 32679, 32647, 32610,
  652.         32568, 32522, 32470, 32413, 32352, 32286, 32214, 32138, 32058,
  653.         31972, 31881, 31786, 31686, 31581, 31471, 31357, 31238, 31114,
  654.         30986, 30853, 30715, 30572, 30425, 30274, 30118, 29957, 29792,
  655.         29622, 29448, 29269, 29086, 28899, 28707, 28511, 28311, 28106,
  656.         27897, 27684, 27467, 27246, 27020, 26791, 26557, 26320, 26078,
  657.         25833, 25583, 25330, 25073, 24812, 24548, 24279, 24008, 23732,
  658.         23453, 23170, 22884, 22595, 22302, 22006, 21706, 21403, 21097,
  659.         20788, 20475, 20160, 19841, 19520, 19195, 18868, 18538, 18205,
  660.         17869, 17531, 17190, 16846, 16500, 16151, 15800, 15447, 15091,
  661.         14733, 14373, 14010, 13646, 13279, 12910, 12540, 12167, 11793,
  662.         11417, 11039, 10660, 10279,  9896,  9512,  9127,  8740,  8351,
  663.         7962,   7571,  7180,  6787,  6393,  5998,  5602,  5205,  4808,
  664.         4410,   4011,  3612,  3212,  2811,  2411,  2009,  1608,  1206,
  665.         804,     402,     0
  666.     };
  668.     if (firstcall) {
  669.         /* Initialization: fill cosine table */
  670.         firstcall = 0;
  672. /* for (i = 0; i < DFTLENGTH; i++ )
  673.      lsp_cos[i] = cos(i*(TWOPI / DFTLENGTH)); */
  674. for (i = 0; i < DFTLENGTH_D4; i++) {
  675.   lsp_cos[i] = cos_tab[i];
  676.   lsp_cos[i+DFTLENGTH_D4] = -cos_tab[DFTLENGTH_D4-i];
  677.   lsp_cos[i+2*DFTLENGTH_D4] = -cos_tab[i];
  678.   lsp_cos[i+3*DFTLENGTH_D4] = cos_tab[DFTLENGTH_D4-i];
  679. }
  681. /* compute default values for w[] */
  682. /* (1./p2) in Q15 */
  683. default_w = divide_s(ONE_Q11,shl_a(p,10));
  684. /* w[0] = (0.5/p2) in Q15 */
  685. /*  default_w0 = shr(default_w,1);   */
  686.     L_temp = _shr2((Longword)default_w,1);
  687.     default_w0 = (Shortword) (0x0000ffffL & L_temp);
  688.     }
  690.     prev_less = 1;
  691.     mag0 = (Longword)0x7fffffff; 
  692. //         L_data_move();mark del
  693.     mag1 = (Longword)0x7fffffff; 
  694. //         L_data_move();mark del
  695.     s_mag0 = (Shortword)0x7fff;  
  696. //        data_move();mark del
  697.     s_mag1 = (Shortword)0x7fff;
  698. //          data_move();mark del
  699.    /* p2 = shr(p,1);       p/2 */
  700.    L_temp = _shr2((Longword)p,1);
  701.    p2 = (Shortword) (0x0000ffffL & L_temp);
  702.     count = 0;      
  703. //    data_move(); mark del
  704.     p_begcos = &lsp_cos[0];
  705.     p_endcos = &lsp_cos[DFTLENGTH-1];
  706.     p_endlsp = &lsp[p2];
  708.     /*  Search all frequencies for minima of Pc(w) */
  709.     for (i = 0; i <= DFTLENGTH_D2; i++ ) {
  710.         p_lsp = p_endlsp;
  711.         p_cos = p_begcos++;
  712. /* mag2 = 0.5 * *(p_lsp--); */
  713. /*  mag2 = L_mult(*(p_lsp--),(Shortword)X05_Q14);      mag2 in Q25 */    
  714.       mag2 = _smpy(*(p_lsp--),(Shortword)X05_Q14);     /* mag2 in Q25 */
  715.         for (j = p2 - 1; j >= 0; j-- ) {
  716.   /*   L_temp1 = L_shr(L_mult(*(p_lsp--),*(p_cos)),1); */
  717. /*             mag2 = L_add(mag2,L_temp1);                 */
  718.             L_temp1 = _sshvr(_smpy(*(p_lsp--),*(p_cos)),1);
  719.             mag2 = _sadd(mag2,L_temp1);
  720. //     complexity -= 4; mark del   /* this can be done with MAC on TMS320C5x */
  721.             p_cos += i;
  722.             if (p_cos > p_endcos)
  723.                 p_cos -= DFTLENGTH;
  724.         }
  725. /* s_mag2 = extract_h(mag2); */
  726. /*         mag2 = L_abs(mag2);   */
  727.      s_mag2 = (Shortword) (0x0000ffffL & (mag2 >> 16));
  728.         mag2 = _abs(mag2);
  730. // L_compare_nonzero();  mark del
  731.         if (mag2 < mag1) {
  732.     prev_less = 1;      
  733. //     data_move();mark del
  734. }
  735.         else {
  736.             if (prev_less) {
  737.       /* check for sign change */
  738. //       logic();    compare_zero();mark del
  739.       if ((s_mag0 ^ s_mag2) < 0) {
  740. /* Minimum frequency found */
  741. /*w[count] = i - 1 + (0.5 *
  742.              (mag0 - mag2) / (mag0 + mag2 - 2*mag1) );
  743.      w[count] *= (2. / DFTLENGTH) ;*/
  744. L_temp1 = _sshvr(L_sub(mag0,mag2),1);
  745. L_temp2 = _ssub(mag0,L_shl(mag1,1));
  746. L_temp2 = _sadd(L_temp2,mag2);
  747. temp1 = L_divider2(L_temp1,L_temp2,0,0);  /* temp1 in Q15 */
  748. /*  temp1 = shr(temp1,9);     Q6 */  
  749.     L_temp = _shr2((Longword)temp1,9);
  750. temp1 = (Shortword) (0x0000ffffL & L_temp);
  751. temp2 = sub(i,1);
  752. temp1 = add(shl_a(temp2,6),temp1);
  753. w[count] = divide_s(temp1,shl_a(DFTLENGTH,5));
  754. //       data_move();mark del
  755. /*  count = add(count,1);  */
  756.     L_temp = _sadd2((Longword)count,(Longword)1);
  757. count = (Shortword) (0x0000ffffL & L_temp);
  758.       }
  759.     }
  760.     prev_less = 0;      
  761. //     data_move();mark del
  762.         }
  763.         mag0 = mag1;      
  764. //        L_data_move();mark del
  765.         mag1 = mag2; 
  766. //             L_data_move();mark del
  767. s_mag0 = s_mag1;     
  768. //  data_move();mark del
  769. s_mag1 = s_mag2;     
  770. //  data_move();mark del
  772.     }
  773.     
  774.     /*  Verify that all roots were found                */
  775.     if (count != p2) {
  776. #if (PRINT)
  777.         printf("ERROR: couldn't find LSP frequencies in frame %d.n",frames);
  778.         printf(" Found were ");
  779.         for (i = 0; i < count; i++) 
  780.             printf(" %10.4f ",(float)w[i]/(1<<15));
  781.         printf("n");
  782. #endif
  783.         /* use default values */
  784.         w[0] = default_w0;
  785.         for (i = 1; i < p2; i++) {
  786.             w[i] = w[i-1] + default_w;
  787.         }
  788.     }
  789. } /* LSP_TO_FREQ */
  790. #undef  DFTLENGTH
  791. #undef  DFTLENGTH_D2
  792. #undef  DFTLENGTH_D4
  795. /* LPC_PRED2REFL
  796.       get refl coeffs from the predictor coeffs
  797.       Input:
  798.          a- the predictor coefficients
  799.          p- the predictor order
  800.       Output:
  801.          k- the reflection coefficients
  802.       Returns:
  803.          energy - energy of residual signal
  804.    Reference:  Markel and Gray, Linear Prediction of Speech
  806.    Q values:
  807.    *a - Q12
  808.    *k - Q15
  809.    energy - Q15
  810. */
  812. Shortword lpc_pred2refl(Shortword *a,Shortword *k,Shortword p)
  813. {
  814.     Shortword *b,*b1, e;
  815.     Shortword energy = ONE_Q15;
  816.     Shortword i,j;
  817.     Shortword shift;
  818.     Shortword temp, n_temp, n_e;
  820. Longword L_temp,L_shift;
  821. Shortword temp1;
  823.     MEM_ALLOC(MALLOC,b,p+1,Shortword);
  824.     MEM_ALLOC(MALLOC,b1,p+1,Shortword);
  826.     /* equate temporary variables (b = a) */
  827.     for(i=1; i <= p; i++) {
  828.       b[i] = a[i]; 
  829. //           data_move();mark del
  830.     }
  832.     /* compute reflection coefficients */
  833.     for(i=p; i > 0; i--)
  834.     {
  835.         k[i] = shl_a(b[i],3);     /* Q12 -> Q15 */
  837.         /* e = sub(ONE_Q15,mult(k[i],k[i])); */
  838. /*         energy = mult(energy,e);          */
  839.          L_temp = _smpy(k[i],k[i]);
  840.  temp1 = (Shortword) (0x0000ffffL & (L_temp >> 16));
  841.  e = sub(ONE_Q15,temp1);
  842.  L_temp = _smpy(energy,e);
  843.  energy = (Shortword) (0x0000ffffL & (L_temp >> 16));
  845.         for(j=1; j < i; j++) {
  846.             b1[j] = b[j];   
  847. //               data_move();mark del
  848. }
  850.         for(j=1; j < i; j++) {
  851.     /* b[j] = (b1[j] - k[i]*b1[i-j])/e; */
  852.    /*  temp = mult(k[i],b1[i-j]);     temp in Q12  */
  853.      L_temp = _smpy(k[i],b1[i-j]);
  854.  temp = (Shortword) (0x0000ffffL & (L_temp >> 16));
  856.     temp = sub(b1[j],temp);
  858.     /* check signs of temp and e before division */
  859.    /*  n_temp = abs_s(temp); */
  860. /*      n_e = abs_s(e);          */
  861. /*      shift = norm_s(n_e);     */
  862.         L_temp = _abs2((Longword)temp);
  863. n_temp = (Shortword) (0x0000ffffL & L_temp);
  864. L_temp = _abs2((Longword)e);
  865. n_e = (Shortword) (0x0000ffffL & L_temp);
  866. L_temp = (Longword)n_e << 16;
  867. L_shift = _norm(L_temp);
  868. if (L_temp == 0){
  869.   L_shift = (Longword)0;
  870. }
  871. shift = (Shortword) (0x0000ffffL & L_shift);
  876.     n_e = shl(n_e,shift);
  877.     if (n_temp > n_e) {
  878.       /* n_temp = shr(n_temp,1); */
  879. /*        shift = add(shift,1);      */
  880.           L_temp = _shr2((Longword)n_temp,1);
  881.   n_temp = (Shortword) (0x0000ffffL & L_temp);
  882.   L_shift = _sadd2((Longword)shift,(Longword)1);
  883.   shift = (Shortword) (0x0000ffffL & L_shift);
  884.     }
  885.     b[j] = shl(divide_s(n_temp,n_e),shift);    /* b[j] in Q12 */
  887.     if ((temp ^ e) < 0)
  888.       b[j] = negate(b[j]);
  889. }
  890.     }
  892.     MEM_FREE(FREE,b1);
  893.     MEM_FREE(FREE,b);
  894.     return(energy);
  895. } /* LPC_PRED2REFL */
  898. /* LPC_LSP2PRED
  899.       get predictor coefficients from the LSPs
  900.    Synopsis: lpc_lsp2pred(w,a,p)
  901.       Input:
  902.          w- the LSPs
  903.          p- the predictor order
  904.       Output:
  905.          a- the predictor coefficients
  906.    Reference:  Kabal and Ramachandran
  907. */
  908. /* Q values:
  909.    w - Q15
  910.    a - Q12
  911.    f - Q25
  912.    c - Q14 */
  914. Shortword lpc_lsp2pred(Shortword *w,Shortword *a,Shortword p)
  915. {
  916.     Shortword i,j,k,p2;
  917.     Shortword c[2];
  918.     Longword **f;
  919.     Longword L_temp;
  921. Longword L_temp1;
  922. Shortword temp;
  924.     /* ensure minimum separation and sort */
  925.     (void)lpc_clmp(w,lsp_delta,p);
  927.     /* p2 = p/2 */
  928.     p2 = shr_a(p,1);
  929. //         data_move();mark del
  931.     MEM_2ALLOC(MALLOC,f,2,p2+1,Longword);
  933.     /* f is Q25 */
  934.     f[0][0] = f[1][0] = (Longword)ONE_Q25;
  935. //        L_data_move();    L_data_move();mark del
  937.     /* -2.0*cos((double)w[1]*M_PI) */
  938.    /*  f[0][1] = L_shr(L_deposit_h(negate(cos_fxp(w[1]))),5); */
  939.      temp = negate(cos_fxp(w[1]));
  940.  L_temp1 = (Longword)temp << 16;
  941.      f[0][1] = _sshvr(L_temp1,5);
  943. //      L_data_move();mark del
  944.    /*  f[1][1] = L_shr(L_deposit_h(negate(cos_fxp(w[2]))),5); */
  945.      temp = negate(cos_fxp(w[2]));
  946.  L_temp1 = (Longword)temp << 16;
  947.      f[1][1] = _sshvr(L_temp1,5);
  948. //      L_data_move();mark del
  950.     k = 3;
  952.     for(i=2; i <= p2; i++)
  953.       { 
  954.         /* c is Q14 */
  955.         /* multiply by 2 is considered as Q15->Q14 */
  956.         c[0] = negate(cos_fxp(w[k++]));
  957. //            data_move();mark del
  958.         c[1] = negate(cos_fxp(w[k++]));
  959. //            data_move();mark del
  960.        
  961.         f[0][i] = f[0][i-2];
  962. //            L_data_move();mark del
  963.         f[1][i] = f[1][i-2]; 
  964. //           L_data_move();mark del
  966.         for(j=i; j >= 2; j--)
  967.           {
  968.     /* f[0][j] += c[0]*f[0][j-1]+f[0][j-2] */
  969.     L_temp = L_mpy_ls(f[0][j-1],c[0]);
  970.     /* L_temp = L_add(L_shl(L_temp,1),f[0][j-2]); */
  971. /*             f[0][j] = L_add(f[0][j],L_temp);       */
  972.         L_temp = _sadd(_sshvl(L_temp,1),f[0][j-2]); 
  973.              f[0][j] = _sadd(f[0][j],L_temp);     
  974.         
  975.  //               L_data_move();mark del
  977.             /* f[1][j] += c[1]*f[1][j-1]+f[1][j-2] */
  978.     L_temp = L_mpy_ls(f[1][j-1],c[1]);
  979.    /*  L_temp = L_add(L_shl(L_temp,1),f[1][j-2]); */
  980. /*             f[1][j] = L_add(f[1][j],L_temp);       */
  981.          L_temp = _sadd(_sshvl(L_temp,1),f[1][j-2]);
  982.             f[1][j] = _sadd(f[1][j],L_temp);
  983. //                L_data_move();mark del
  984.           }
  986.        /*  f[0][1] = L_add(f[0][1],L_shl(L_mpy_ls(f[0][0],c[0]),1)); */
  987. /*         f[1][1] = L_add(f[1][1],L_shl(L_mpy_ls(f[1][0],c[1]),1)); */
  988.          f[0][1] = _sadd(f[0][1],_sshl(L_mpy_ls(f[0][0],c[0]),1));
  989.          f[1][1] = _sadd(f[1][1],_sshl(L_mpy_ls(f[1][0],c[1]),1));
  990.       }
  992.     for(i=p2; i > 0; i--)
  993.       {
  994.         /* short f (f is a Q14) */
  995.        /*  f[0][i] = L_add(f[0][i],f[0][i-1]); */
  996. /*         f[1][i] = L_sub(f[1][i],f[1][i-1]); */
  997.           f[0][i] = _sadd(f[0][i],f[0][i-1]);
  998.           f[1][i] = _ssub(f[1][i],f[1][i-1]);
  1000.         /* a is Q12 */
  1001.         /* a[i] = 0.50*(f[0][i]+f[1][i]) */
  1002. /* a[p+1-i] = 0.50*(f[0][i]-f[1][i]) */
  1003. /* Q25 -> Q27 -> Q12 */
  1004.        /*  a[i] = extract_h(L_shl(L_add(f[0][i],f[1][i]),2));     */
  1005. /*         a[p+1-i] = extract_h(L_shl(L_sub(f[0][i],f[1][i]),2)); */
  1006.         L_temp1 = _sshl(_sadd(f[0][i],f[1][i]),2);
  1007.         a[i] = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  1008. L_temp1 = _sshl(_ssub(f[0][i],f[1][i]),2);
  1009. a[p+1-i] = (Shortword) (0x0000ffffL & (L_temp1 >> 16));
  1011.     }
  1013.     MEM_2FREE(FREE,f);
  1014.     return((Shortword)0);
  1015. } /* LPC_LSP2PRED */
  1018. /*
  1019.     Name: lpc_syn- LPC synthesis filter.
  1020.     Aliases: lpc_synthesis
  1021.     Description:
  1022.         LPC all-pole synthesis filter
  1024.         for j = 0 to n-1
  1025.             y[j] = x[j] - sum(k=1 to p) y[j-k] a[k]
  1027.     Inputs:
  1028.         x- input vector (n samples, x[0..n-1])
  1029.         a- lpc vector (order p, a[1..p])
  1030.         p- order of lpc filter
  1031.         n- number of elements in vector which is to be filtered
  1032.         y[-p..-1]- filter memory (past outputs)
  1033.     Outputs:
  1034.         y- output vector (n samples, y[0..n-1])
  1035.     Returns: NULL
  1036.     Includes:
  1037.         spbstd.h
  1038.         lpc.h
  1040.     Systems and Info. Science Lab
  1041.     Copyright (c) 1995 by Texas Instruments, Inc.  All rights reserved.
  1043. */
  1045. Shortword lpc_syn(Shortword *x,Shortword *y,Shortword *a,Shortword p,
  1046.   Shortword n)
  1047. {
  1048.     Shortword i,j;
  1049.     Longword accum;
  1051. Longword L_temp,L_Prod;
  1053.     for(j=0; j < n; j++) {
  1054.         /* accum = L_shr(L_deposit_h(x[j]),3); */
  1055. L_temp = (Longword)(x[j]) << 16;
  1056. accum = _sshvr(L_temp,3);
  1057. for(i=p; i > 0; i--)
  1058.     accum = L_msu(accum,y[j-i],a[i]);
  1059. /* Round off output */
  1060. /* accum = L_shl(accum,3); */
  1061. /* y[j] = round(accum); */
  1062.   accum = _sshl(accum,3);
  1063.   L_Prod = _sadd(accum, 0x00008000L); /* round MSP */
  1064.       y[j] = (Shortword)(0x0000ffffL & (L_Prod >> 16));
  1066.     }
  1067.     return((Shortword)0);
  1068. }