1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
|
/*
libaudiomask - hybrid simultaneous audio masking threshold evaluation library
Copyright (C) 2000-2010 Dr Matthew Raphael Flax
This program is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
/*
*
* This example shows how to use Dr M.R. Flax's (2000) hybrid
* simultaneous audio masking class to find the masking threshold of a time domain signal.
*
* The compilation of this file is demonstrated in Makefile.
* Run this file : ./AudioMaskerExample
* View the results of this file using www.octave.org by running the script view.m
* - simply type view once octave has started and you are in the suitable directory.
*
* The input audio should be stored in the file INPUTFILENAME in text format - each sample seperated by a * white space.
*
* ========================= HOWTO ===============================
* \code
* // First find the masking threshold
* AudioMasker masker(sampleFreq, count); // Create the audio masker class using fs=sampleFreq and count filters
* masker.excite(input, sampleCount); // find the mask for the array of input data which has sampleCount time samples.
*
* // Now do something with the masking threshold ...
*
* // The frequency domain mask is now located here
* for (int j=0; j<count;j++)
* masker.mask[j]; // This is the mask at each of the count frequencies of interest
*
* // A more sophisticated example - find the threshold for each Fourier bin
* double fact=(double)sampleFreq/((double)sampleCount-1.0); // convert from an index to the equivalent * Fourier bin frequency
* for (int j=0; j<halfSampleCount;j++){
* cout<<"finding for freq "<<j*fact<<'\t'; // The frequency we are inspecting
* double threshold=masker.findThreshold(j*fact); // The masking threshold
* 20*log10(threshold); // The threshold in decibels (dB)
* }
*
* \endcode
* // The following example calculates the spectral power and the masking threshold for each Fourier bin of interest ...
*/
#include <math.h>
#include "AudioMasker.H"
#define INPUTFILENAME "audio.44100.txt" // input file - text written samples seperated by white spaces
#define POWFILENAME "fa.pow" // The power spectrum of the input signal
#define THRESHFILENAME "thresh.dat" // The masking threshold at each frequency of the input signal
#define TMASKFILENAME "fa.t.mask" // The masking threshold for each filter CF
#define EXCITEFILENAME "fa.excite" // The excitation
#include <fstream>
int main(void){
// Setup many variables
// The number of time domain samples
int sampleCount=1024, halfSampleCount=(int)rint((double)sampleCount/2.0);
// The filter bank count and sample frequency
int count=50, sampleFreq=44100;
// The number of time domain samples to skip before the sample of interest.
int skip=8192-sampleCount-1;
// The input array to hold the input time data
double input[sampleCount];
// open the input file
ifstream inputF(INPUTFILENAME);
// Skip the first 2*'skip' samples.
int temp;
for (int i=0; i<skip;i++)
inputF >> temp >> input[0];
// load sampleCount samples as the input to the algorithm
for (int i=0; i<sampleCount;i++)
inputF >> input[i];
inputF.close();
ofstream outputCF("cf.dat"); // central freq. output file
ofstream outputT(TMASKFILENAME);
ofstream outputP(POWFILENAME); // Input data Fourier power file
// Get our masking function (class) ...
AudioMasker masker(sampleFreq, count);
//AudioMasker masker; // Can also be called like so with default filter banks and sampleFrequency
masker.excite(input, sampleCount); // find the mask
for (int j=0; j<count;j++){ // Output the central freq to file
outputCF <<masker.pfb->cf[j]*((double)sampleCount/(double)sampleFreq)<<'\t';
outputT << 20*log10(masker.mask[j])<<'\t'; // output the mask for each filter CF to file
}
outputCF<<endl;
outputT<<endl;
realFFTData fftData(sampleCount); // Find the fourier transform for the output of the power to file
realFFT fft(&fftData); // init the Fourier transform
for (int j=0; j<sampleCount;j++) // load the time domain input to the Fourier transform
fftData.in[j]=input[j];
fft.fwdTransform(); // Perform the transformation of the time domain to the frequency domain.
fftData.compPowerSpec(); // Find the power spectrum
for (int j=0; j<sampleCount/2;j++) // Output the power spectrum to file
outputP<<20*log10(sqrt(fftData.power_spectrum[j]))<<'\t';
outputP<<endl;
outputCF.close();
outputT.close();
outputP.close();
ofstream outputF(THRESHFILENAME); // Find and output the masking threshold for each Fourier bin in the power spectrum
double fact=(double)sampleFreq/((double)sampleCount-1.0); // convert from an index to the equivalent Fourier bin frequency
for (int j=0; j<halfSampleCount;j++){
// cout<<"finding for freq "<<j*fact<<'\t';
outputF<<20*log10(masker.findThreshold(j*fact))<<'\t';
}
outputF<<endl;
outputF.close();
}
|
