1 /*-------------- Telecommunications & Signal Processing Lab ---------------
2                              McGill University
3 
4 Routine:
5   struct MOVB_MDiff PQmovModDiffB (const double *M[2], const double ERavg[],
6                                    const struct Par_MDiff *MDiff)
7   void PQavgModDiffB (int Ndel, const struct MOVC_MDiff *MDiff, int Nchan,
8                       int Np, double *WinModDiff1B, double *AvgModDiff1B,
9 	              double *AvgModDiff2B)
10 
11 Purpose:
12   Calculate modulation difference MOV precursors (PEAQ Basic version)
13   Average modulation differences (PEAQ Basic version)
14 
15 Description:
16   PQmovModDiffB:
17   This routine calculates the modulation difference MOV precursors.
18 
19   PQavgModDiffB:
20   This routine averages the modulation difference values to give the
21   final MOV values.
22 
23 Parameters:
24   PQmovModDiffB:
25   <-  struct MOVB_MDiff PQmovModB
26       Output modulation differene MOV precursors
27    -> const double *M[2]
28       Modulation patterns
29    -> const double ERavg[]
30       Modified excitation patterns
31    -> const struct Par_MDiff *MDiff
32       Modulation difference processing parameters
33 
34   PQavgModDiffB:
35    -> int Ndel
36       Number of PEAQ frames to skip at the start
37    -> const struct MOVC_MDiff *MDiff
38       Modulation difference values
39    -> int Nchan
40       Number of channels
41    -> int Np
42       Number of PEAQ frames
43   <-  double *WinModDiff1B
44       Model output variable
45   <-  double *AvgModDiff1B
46       Model output variable
47   <-  double *AvgModDiff2B
48       Model output variable
49 
50 Author / revision:
51   P. Kabal  Copyright (C) 2003
52   $Revision: 1.5 $  $Date: 2003/05/13 01:15:12 $
53 
54 -------------------------------------------------------------------------*/
55 
56 #include <math.h>
57 
58 #include "PQevalAudio.h"
59 
60 static double
61 PQ_WinAvg (int L, const double x[], int N);
62 static double
63 PQ_WtAvg (const double x[], const double W[], int N);
64 
65 
66 struct MOVB_MDiff
PQmovModDiffB(const double * M[2],const double ERavg[],const struct Par_MDiff * MDiff)67 PQmovModDiffB (const double *M[2], const double ERavg[],
68 	       const struct Par_MDiff *MDiff)
69 
70 {
71   int m;
72   double s1B, s2B, Wt, num1B, num2B, MD1B, MD2B;
73   struct MOVB_MDiff MDiffV;
74 
75   const int Nc = MDiff->Nc;
76   const double *Ete = MDiff->Ete;
77 
78   /* Parameters */
79   static const double negWt2B = 0.1;
80   static const double offset1B = 1.0;
81   static const double offset2B = 0.01;
82   static const double levWt = 100;
83 
84 /* Mt1B - Modulation Difference 1B
85    Mt2B - Modulation Difference 2B
86    Wt   - Weighting value
87 */
88   s1B = 0;
89   s2B = 0;
90   Wt = 0;
91   for (m = 0; m < Nc; ++m) {
92     if (M[0][m] > M[1][m]) {
93       num1B = M[0][m] - M[1][m];
94       num2B = negWt2B * num1B;
95     }
96     else {
97       num1B = M[1][m] - M[0][m];
98       num2B = num1B;
99     }
100     MD1B = num1B / (offset1B + M[0][m]);
101     MD2B = num2B / (offset2B + M[0][m]);
102     s1B += MD1B;
103     s2B += MD2B;
104     Wt += ERavg[m] / (ERavg[m] + levWt * Ete[m]);
105   }
106 
107   MDiffV.Mt1B = (100. / Nc) * s1B;
108   MDiffV.Mt2B = (100. / Nc) * s2B;
109   MDiffV.Wt = Wt;
110 
111   return MDiffV;
112 }
113 
114 void
PQavgModDiffB(int Ndel,const struct MOVC_MDiff * MDiff,int Nchan,int Np,double * WinModDiff1B,double * AvgModDiff1B,double * AvgModDiff2B)115 PQavgModDiffB (int Ndel, const struct MOVC_MDiff *MDiff, int Nchan, int Np,
116 	       double *WinModDiff1B, double *AvgModDiff1B,
117 	       double *AvgModDiff2B)
118 
119 {
120   int j, L;
121   double s;
122 
123   const double tavg = 0.1;
124   const double Fss = PQ_FS / PQ_CB_ADV;
125 
126   /* Windowed average - delayed average */
127   L = (int) floor (tavg * Fss);		/* 100 ms sliding window length */
128   s = 0;
129   for (j = 0; j < Nchan; ++j)
130     s += PQ_WinAvg (L, &MDiff->Mt1B[j][Ndel], Np-Ndel);
131   *WinModDiff1B = s / Nchan;
132 
133   /* Weighted linear average - delayed average */
134   s = 0;
135   for (j = 0; j < Nchan; ++j)
136     s += PQ_WtAvg (&MDiff->Mt1B[j][Ndel], &MDiff->Wt[j][Ndel], Np-Ndel);
137   *AvgModDiff1B = s / Nchan;
138 
139   /* Weighted linear average - delayed average */
140   s = 0;
141   for (j = 0; j < Nchan; ++j)
142     s += PQ_WtAvg (&MDiff->Mt2B[j][Ndel], &MDiff->Wt[j][Ndel], Np-Ndel);
143   *AvgModDiff2B = s / Nchan;
144 
145   return;
146 }
147 
148 /* Sliding window (L samples) average */
149 
150 static double
PQ_WinAvg(int L,const double x[],int N)151 PQ_WinAvg (int L, const double x[], int N)
152 {
153   int i, m;
154   double s, t;
155 
156   s = 0;
157   for (i = L-1; i < N; ++i) {
158     t = 0;
159     for (m = 0; m < L; ++m)
160       t += sqrt (x[i-m]);
161     s += pow (t / L, 4);
162   }
163 
164   if (N >= L)
165     s = sqrt (s / (N - L + 1));
166 
167   return s;
168 }
169 
170 /* Weighted average */
171 
172 static double
PQ_WtAvg(const double x[],const double W[],int N)173 PQ_WtAvg (const double x[], const double W[], int N)
174 
175 {
176   int i;
177   double s;
178   double sW;
179 
180   s = 0;
181   sW = 0;
182   for (i = 0; i < N; ++i) {
183     s += W[i] * x[i];
184     sW += W[i];
185   }
186 
187   if (N > 0)
188     s = s / sW;
189 
190   return s;
191 }
192