resample-test.cc
Go to the documentation of this file.
1 // feat/resample-test.cc
2 
3 // Copyright 2013 Pegah Ghahremani
4 // 2014 IMSL, PKU-HKUST (author: Wei Shi)
5 // 2014 Yanqing Sun, Junjie Wang
6 // 2014 Johns Hopkins University (author: Daniel Povey)
7 
8 // See ../../COPYING for clarification regarding multiple authors
9 //
10 // Licensed under the Apache License, Version 2.0 (the "License");
11 // you may not use this file except in compliance with the License.
12 // You may obtain a copy of the License at
13 //
14 // http://www.apache.org/licenses/LICENSE-2.0
15 //
16 // THIS CODE IS PROVIDED *AS IS* BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
17 // KIND, EITHER EXPRESS OR IMPLIED, INCLUDING WITHOUT LIMITATION ANY IMPLIED
18 // WARRANTIES OR CONDITIONS OF TITLE, FITNESS FOR A PARTICULAR PURPOSE,
19 // MERCHANTABLITY OR NON-INFRINGEMENT.
20 // See the Apache 2 License for the specific language governing permissions and
21 // limitations under the License.
22 
23 
24 #include "feat/resample.h"
25 
26 using namespace kaldi;
27 
28 class TestFunction {
29  public:
30  explicit TestFunction(double frequency):
31  frequency_(frequency),
32  sin_magnitude_(RandGauss()),
33  cos_magnitude_(RandGauss()) { }
34 
35  double operator() (double t) const {
36  double omega_t = t * M_2PI * frequency_;
37  return sin_magnitude_ * sin(omega_t)
38  + cos_magnitude_ * cos(omega_t);
39  }
40  private:
41  double frequency_;
44 };
45 
46 
48  BaseFloat samp_freq = 1000.0 * (1.0 + RandUniform());
49  int32 num_samp = 256 + static_cast<int32>((RandUniform() * 256));
50 
51  BaseFloat time_interval = num_samp / samp_freq;
52 
53  // Choose a lowpass frequency that's lower than 95% of the Nyquist.
54  BaseFloat lowpass_freq = samp_freq * 0.95 * 0.5 / (1.0 + RandUniform());
55 
56  // Number of zeros of the sinc function that the window extends out to.
57  int32 num_zeros = 3 + rand() % 10;
58 
59  // Resample the signal at arbitrary points within that time interval.
60  int32 num_resamp = 50 + rand() % 100; // Resample at around 100 points,
61  // anywhere in the signal.
62 
63 
64  Vector<BaseFloat> resample_points(num_resamp);
65  for (int32 i = 0; i < num_resamp; i++) {
66  // the if-statement is to make some of the resample_points
67  // exactly coincide with the original points, to activate
68  // a certain code path.
69  if (rand() % 2 == 0)
70  resample_points(i) = (rand() % num_samp) / samp_freq;
71  else
72  resample_points(i) = RandUniform() * time_interval;
73  }
74 
75 
76 
77  BaseFloat window_width = num_zeros / (2.0 * lowpass_freq);
78  // the resampling should be quite accurate if we are further
79  // than filter_width away from the edges.
80  BaseFloat min_t = 0.0 + window_width,
81  max_t = time_interval - (1.0 / samp_freq) - window_width;
82 
83  // window_freq gives us a rough idea of the frequency spread
84  // that the windowing function gives us; we want the test frequency
85  // to be lower than the lowpass frequency by at least this much.
86  // (note: the real width of the window from side to side
87  // is 2.0 * window_width)
88  BaseFloat window_freq = 1.0 / (2.0 * window_width),
89  freq_margin = 2.0 * window_freq;
90 
91  // Choose a test-signal frequency that's lower than
92  // lowpass_freq - freq_margin.
93  BaseFloat test_signal_freq =
94  (lowpass_freq - freq_margin) * (1.0 / (1.0 + RandUniform()));
95 
96  KALDI_ASSERT(test_signal_freq > 0.0);
97 
98  ArbitraryResample resampler(num_samp, samp_freq, lowpass_freq,
99  resample_points, num_zeros);
100 
101 
102  TestFunction test_func(test_signal_freq);
103 
104  // test with a one-row matrix equal to the test signal.
105  Matrix<BaseFloat> sample_values(1, num_samp);
106  for (int32 i = 0; i < num_samp; i++) {
107  BaseFloat t = i / samp_freq;
108  sample_values(0, i) = test_func(t);
109  }
110  Matrix<BaseFloat> resampled_values(1, num_resamp);
111 
112 
113  if (rand() % 2 == 0) {
114  resampler.Resample(sample_values,
115  &resampled_values);
116  } else {
117  SubVector<BaseFloat> out(resampled_values, 0);
118  resampler.Resample(sample_values.Row(0),
119  &out);
120  }
121 
122 
123  for (int32 i = 0; i < num_resamp; i++) {
124  BaseFloat t = resample_points(i),
125  x1 = test_func(t),
126  x2 = resampled_values(0, i),
127  error = fabs(x1 - x2);
128  if (i % 10 == 0) {
129  KALDI_VLOG(1) << "Error is " << error << ", t = " << t
130  << ", samp_freq = " << samp_freq << ", lowpass_freq = "
131  << lowpass_freq << ", test_freq = " << test_signal_freq
132  << ", num-zeros is " << num_zeros;
133  }
134  if (t > min_t && t < max_t) {
135  if (num_zeros == 3) {
136  KALDI_ASSERT(error < 0.1);
137  } else {
138  KALDI_ASSERT(error < 0.025);
139  }
140  } else {
141  KALDI_VLOG(1) << "[not checking since out of bounds]";
142  }
143  }
144 }
145 
146 
148  // this test makes sure that LinearResample gives identical results to
149  // ArbitraryResample when set up the same way, even if the signal is broken up
150  // into many pieces.
151 
152  int32 samp_freq = 1000.0 * (1.0 + RandUniform()),
153  resamp_freq = 1000.0 * (1.0 + RandUniform());
154  // note: these are both integers!
155  int32 num_samp = 256 + static_cast<int32>((RandUniform() * 256));
156 
157  BaseFloat time_interval = num_samp / static_cast<BaseFloat>(samp_freq);
158 
159  // Choose a lowpass frequency that's lower than 95% of the Nyquist of both
160  // of the frequencies..
161  BaseFloat lowpass_freq =
162  std::min(samp_freq, resamp_freq) * 0.95 * 0.5 / (1.0 + RandUniform());
163 
164  // Number of zeros of the sinc function that the window extends out to.
165  int32 num_zeros = 3 + rand() % 10;
166 
167  // compute the number of "resample" points.
168  int32 num_resamp = ceil(time_interval * resamp_freq);
169 
170  Vector<BaseFloat> resample_points(num_resamp);
171  for (int32 i = 0; i < num_resamp; i++)
172  resample_points(i) = i / static_cast<BaseFloat>(resamp_freq);
173 
174 
175  Vector<BaseFloat> test_signal(num_samp);
176  test_signal.SetRandn();
177 
178  ArbitraryResample resampler(num_samp, samp_freq, lowpass_freq,
179  resample_points, num_zeros);
180 
181 
182  // test with a one-row matrix equal to the test signal.
183  Matrix<BaseFloat> sample_values(1, num_samp);
184  sample_values.Row(0).CopyFromVec(test_signal);
185 
186  Matrix<BaseFloat> resampled_values(1, num_resamp);
187 
188  resampler.Resample(sample_values,
189  &resampled_values);
190 
191  LinearResample linear_resampler(samp_freq, resamp_freq,
192  lowpass_freq, num_zeros);
193 
194  Vector<BaseFloat> resampled_vec;
195 
196  linear_resampler.Resample(test_signal, true, &resampled_vec);
197 
198  if (!ApproxEqual(resampled_values.Row(0), resampled_vec)) {
199  KALDI_LOG << "ArbitraryResample: " << resampled_values.Row(0);
200  KALDI_LOG << "LinearResample: " << resampled_vec;
201  KALDI_ERR << "Signals differ.";
202  }
203 
204  // Check it gives the same results when the input is broken up into pieces.
205  Vector<BaseFloat> resampled_vec2;
206  int32 input_dim_seen = 0;
207  while (input_dim_seen < test_signal.Dim()) {
208  int32 dim_remaining = test_signal.Dim() - input_dim_seen;
209  int32 piece_size = rand() % std::min(dim_remaining + 1, 10);
210  KALDI_VLOG(1) << "Piece size = " << piece_size;
211  SubVector<BaseFloat> in_piece(test_signal, input_dim_seen, piece_size);
212  Vector<BaseFloat> out_piece;
213  bool flush = (piece_size == dim_remaining);
214  linear_resampler.Resample(in_piece, flush, &out_piece);
215  int32 old_output_dim = resampled_vec2.Dim();
216  resampled_vec2.Resize(old_output_dim + out_piece.Dim(), kCopyData);
217  resampled_vec2.Range(old_output_dim, out_piece.Dim())
218  .CopyFromVec(out_piece);
219  input_dim_seen += piece_size;
220  }
221 
222  if (!ApproxEqual(resampled_values.Row(0), resampled_vec2)) {
223  KALDI_LOG << "ArbitraryResample: " << resampled_values.Row(0);
224  KALDI_LOG << "LinearResample[broken-up]: " << resampled_vec2;
225  KALDI_ERR << "Signals differ.";
226  }
227 }
228 
230  int32 num_samp = 150 + rand() % 100;
231  BaseFloat samp_freq = 1000, resamp_freq = 4000;
232 
233  int32 num_zeros = 10; // fairly accurate.
234  Vector<BaseFloat> signal_orig(num_samp);
235  signal_orig.SetRandn();
236 
237  Vector<BaseFloat> signal(num_samp);
238  { // make sure signal is sufficiently low pass, i.e. that we have enough
239  // headroom before the Nyquist.
240  LinearResample linear_resampler_filter(samp_freq, samp_freq,
241  0.8 * samp_freq / 2.0, num_zeros);
242  linear_resampler_filter.Resample(signal_orig, true, &signal);
243  }
244 
245 
246  Vector<BaseFloat> signal_upsampled;
247 
248  LinearResample linear_resampler(samp_freq, resamp_freq,
249  samp_freq / 2.0, num_zeros);
250 
251  linear_resampler.Resample(signal, true, &signal_upsampled);
252 
253  // resample back to the original frequency.
254  LinearResample linear_resampler2(resamp_freq, samp_freq,
255  samp_freq / 2.0, num_zeros);
256 
257 
258  Vector<BaseFloat> signal_downsampled;
259  linear_resampler2.Resample(signal_upsampled, true, &signal_downsampled);
260 
261 
262  int32 samp_discard = 30; // Discard 20 samples for edge effects.
263  SubVector<BaseFloat> signal_middle(signal, samp_discard,
264  signal.Dim() - (2 * samp_discard));
265 
266  SubVector<BaseFloat> signal2_middle(signal_downsampled, samp_discard,
267  signal.Dim() - (2 * samp_discard));
268 
269  BaseFloat self1 = VecVec(signal_middle, signal_middle),
270  self2 = VecVec(signal2_middle, signal2_middle),
271  cross = VecVec(signal_middle, signal2_middle);
272  KALDI_LOG << "Self1 = " << self1 << ", self2 = " << self2
273  << ", cross = " << cross;
274  AssertEqual(self1, self2, 0.001);
275  AssertEqual(self1, cross, 0.001);
276 }
277 
278 int main() {
279  try {
280  for (int32 x = 0; x < 50; x++)
282  for (int32 x = 0; x < 50; x++)
284  for (int32 x = 0; x < 50; x++)
286 
287  KALDI_LOG << "Tests succeeded.\n";
288  return 0;
289  } catch(const std::exception &e) {
290  KALDI_ERR << e.what();
291  return 1;
292  }
293 }
double frequency_
This code computes Goodness of Pronunciation (GOP) and extracts phone-level pronunciation feature for...
Definition: chain.dox:20
void Resample(const VectorBase< BaseFloat > &input, bool flush, Vector< BaseFloat > *output)
This function does the resampling.
Definition: resample.cc:152
float RandUniform(struct RandomState *state=NULL)
Returns a random number strictly between 0 and 1.
Definition: kaldi-math.h:151
double sin_magnitude_
double cos_magnitude_
void UnitTestArbitraryResample()
float RandGauss(struct RandomState *state=NULL)
Definition: kaldi-math.h:155
void UnitTestLinearResample2()
kaldi::int32 int32
void Resize(MatrixIndexT length, MatrixResizeType resize_type=kSetZero)
Set vector to a specified size (can be zero).
float BaseFloat
Definition: kaldi-types.h:29
const SubVector< Real > Row(MatrixIndexT i) const
Return specific row of matrix [const].
Definition: kaldi-matrix.h:188
Class ArbitraryResample allows you to resample a signal (assumed zero outside the sample region...
Definition: resample.h:95
void Resample(const MatrixBase< BaseFloat > &input, MatrixBase< BaseFloat > *output) const
This function does the resampling.
Definition: resample.cc:280
#define KALDI_ERR
Definition: kaldi-error.h:147
MatrixIndexT Dim() const
Returns the dimension of the vector.
Definition: kaldi-vector.h:64
void SetRandn()
Set vector to random normally-distributed noise.
int main()
TestFunction(double frequency)
A class representing a vector.
Definition: kaldi-vector.h:406
LinearResample is a special case of ArbitraryResample, where we want to resample a signal at linearly...
Definition: resample.h:147
#define KALDI_ASSERT(cond)
Definition: kaldi-error.h:185
static void AssertEqual(float a, float b, float relative_tolerance=0.001)
assert abs(a - b) <= relative_tolerance * (abs(a)+abs(b))
Definition: kaldi-math.h:276
#define M_2PI
Definition: kaldi-math.h:52
#define KALDI_VLOG(v)
Definition: kaldi-error.h:156
void UnitTestLinearResample()
#define KALDI_LOG
Definition: kaldi-error.h:153
Real VecVec(const VectorBase< Real > &a, const VectorBase< Real > &b)
Returns dot product between v1 and v2.
Definition: kaldi-vector.cc:37
Represents a non-allocating general vector which can be defined as a sub-vector of higher-level vecto...
Definition: kaldi-vector.h:501
static bool ApproxEqual(float a, float b, float relative_tolerance=0.001)
return abs(a - b) <= relative_tolerance * (abs(a)+abs(b)).
Definition: kaldi-math.h:265
SubVector< Real > Range(const MatrixIndexT o, const MatrixIndexT l)
Returns a sub-vector of a vector (a range of elements).
Definition: kaldi-vector.h:94