1 /*-------------- Telecommunications & Signal Processing Lab ---------------
2                              McGill University
3 
4 Routine:
5   double PQmovNLoudB (const double *M[2], const double *EP[2],
6                       const struct Par_NLoud *NLoud);
7   double PQavgNLoud (int Ndel, const struct MOVC_Mod *Mod, int Nchan, int Np)
8 
9 Purpose:
10   Calculate the noise loudness MOV precursor (PEAQ Basic version)
11   Calculate the average noise loudness (PEAQ Basic Version)
12 
13 Description:
14   PQmovNLoudB:
15   This routine calculates the modulation difference MOV precursors.
16 
17   PQavgNLoudB:
18   This rouitne averages the modulation difference values.
19 
20 Parameters:
21   PQmovNLoudB:
22   <-  double PQmovNLoudB
23       Output noise loudness MOV precursor
24    -> const double *M[2]
25       Modulation patterns
26    -> const double *EP[2]
27       Excitation patterns
28    -> const struct Par_NLoud *NLoud
29       Modulation difference processing parameters
30 
31   PQavgNLoudB:
32   <-  double PQavgNLoudB
33       Model output variable
34    -> int Ndel
35       Number of PEAQ frames to skip at the start
36    -> const struct MOVC_NLoud *NLoud
37       Noise loudness values
38    -> int Nchan
39       Number of channels
40    -> int Np
41       Number of PEAQ frames
42 
43 Author / revision:
44   P. Kabal  Copyright (C) 2003
45   $Revision: 1.4 $  $Date: 2003/05/13 01:15:12 $
46 
47 -------------------------------------------------------------------------*/
48 
49 #include <math.h>
50 
51 #include "PQevalAudio.h"
52 
53 #define MAXV(a, b)	(((a) > (b)) ? (a) : (b))
54 #define SQRV(x)		((x) * (x))
55 
56 static double
57 PQ_RMSAvg (const double x[], int N);
58 
59 
60 double
PQmovNLoudB(const double * M[2],const double * EP[2],const struct Par_NLoud * NLoud)61 PQmovNLoudB (const double *M[2], const double *EP[2],
62 	     const struct Par_NLoud *NLoud)
63 
64 {
65   int m;
66   double s, sref, stest, beta, a, b, NL;
67 
68   const int Nc = NLoud->Nc;
69   const double *Et = NLoud->Et;
70 
71   static const double alpha = 1.5;
72   static const double TF0 = 0.15;
73   static const double S0 = 0.5;
74   static const double NLmin = 0;
75   static const double e = 0.23;
76 
77   s = 0;
78   for (m = 0; m < Nc; ++m) {
79     sref  = TF0 * M[0][m] + S0;
80     stest = TF0 * M[1][m] + S0;
81     beta = exp (-alpha * (EP[1][m] - EP[0][m]) / EP[0][m]);
82     a = MAXV (0, stest * EP[1][m] - sref * EP[0][m]);
83     b = Et[m] + sref * EP[0][m] * beta;
84     s += pow (Et[m] / stest, e) * (pow (1 + a / b, e) - 1);
85   }
86 
87   NL = (24. / Nc) * s;
88   if (NL < NLmin)
89     NL = 0;
90 
91   return NL;
92 }
93 
94 /* Average noise loudness values */
95 
96 double
PQavgNLoud(int Ndel,const struct MOVC_NLoud * Nloud,int Nchan,int Np)97 PQavgNLoud (int Ndel, const struct MOVC_NLoud *Nloud, int Nchan, int Np)
98 
99 {
100   int j;
101   double s, RmsNoiseLoudB;
102 
103   /* RMS average - delayed average and loudness threshold */
104   s = 0;
105   for (j = 0; j < Nchan; ++j)
106     s += PQ_RMSAvg (&Nloud->NL[j][Ndel], Np-Ndel);
107   RmsNoiseLoudB = s / Nchan;
108 
109   return RmsNoiseLoudB;
110 }
111 
112 /* Square root of average of squared values */
113 
114 static double
PQ_RMSAvg(const double x[],int N)115 PQ_RMSAvg (const double x[], int N)
116 
117 {
118   int i;
119   double s;
120 
121   s = 0;
122   for (i = 0; i < N; ++i)
123     s += SQRV(x[i]);
124 
125   if (N > 0)
126     s = sqrt (s / N);
127 
128   return s;
129 }
130