Merge beamformer
R=andrew@webrtc.org
Review URL: https://webrtc-codereview.appspot.com/34529004
git-svn-id: http://webrtc.googlecode.com/svn/trunk@7958 4adac7df-926f-26a2-2b94-8c16560cd09d
diff --git a/webrtc/modules/audio_processing/BUILD.gn b/webrtc/modules/audio_processing/BUILD.gn
index 01a4d7b..da121f3 100644
--- a/webrtc/modules/audio_processing/BUILD.gn
+++ b/webrtc/modules/audio_processing/BUILD.gn
@@ -72,6 +72,10 @@
"audio_buffer.h",
"audio_processing_impl.cc",
"audio_processing_impl.h",
+ "beamformer/complex_matrix.h",
+ "beamformer/covariance_matrix_generator.cc",
+ "beamformer/covariance_matrix_generator.h",
+ "beamformer/matrix.h",
"channel_buffer.cc",
"channel_buffer.h",
"common.h",
@@ -167,6 +171,13 @@
]
}
+ if (rtc_use_openmax_dl) {
+ sources += [
+ "beamformer/beamformer.cc",
+ "beamformer/beamformer.h",
+ ]
+ }
+
if (cpu_arch == "x86" || cpu_arch == "x64") {
deps += [ ":audio_processing_sse2" ]
}
diff --git a/webrtc/modules/audio_processing/audio_processing.gypi b/webrtc/modules/audio_processing/audio_processing.gypi
index 70c549c..a2c26a3 100644
--- a/webrtc/modules/audio_processing/audio_processing.gypi
+++ b/webrtc/modules/audio_processing/audio_processing.gypi
@@ -8,12 +8,6 @@
{
'variables': {
- 'audio_processing_dependencies': [
- '<(webrtc_root)/base/base.gyp:rtc_base_approved',
- '<(webrtc_root)/common_audio/common_audio.gyp:common_audio',
- '<(webrtc_root)/modules/modules.gyp:iSAC',
- '<(webrtc_root)/system_wrappers/source/system_wrappers.gyp:system_wrappers',
- ],
'shared_generated_dir': '<(SHARED_INTERMEDIATE_DIR)/audio_processing/asm_offsets',
},
'targets': [
@@ -31,7 +25,10 @@
'aec_untrusted_delay_for_testing%': 0,
},
'dependencies': [
- '<@(audio_processing_dependencies)',
+ '<(webrtc_root)/base/base.gyp:rtc_base_approved',
+ '<(webrtc_root)/common_audio/common_audio.gyp:common_audio',
+ '<(webrtc_root)/modules/modules.gyp:iSAC',
+ '<(webrtc_root)/system_wrappers/source/system_wrappers.gyp:system_wrappers',
],
'sources': [
'aec/aec_core.c',
@@ -84,6 +81,10 @@
'audio_buffer.h',
'audio_processing_impl.cc',
'audio_processing_impl.h',
+ 'beamformer/complex_matrix.h',
+ 'beamformer/covariance_matrix_generator.cc',
+ 'beamformer/covariance_matrix_generator.h',
+ 'beamformer/matrix.h',
'channel_buffer.cc',
'channel_buffer.h',
'common.h',
@@ -178,6 +179,12 @@
'ns/windows_private.h',
],
}],
+ ['rtc_use_openmax_dl==1', {
+ 'sources': [
+ 'beamformer/beamformer.cc',
+ 'beamformer/beamformer.h',
+ ],
+ }],
['target_arch=="ia32" or target_arch=="x64"', {
'dependencies': ['audio_processing_sse2',],
}],
diff --git a/webrtc/modules/audio_processing/audio_processing_tests.gypi b/webrtc/modules/audio_processing/audio_processing_tests.gypi
index 7c36be7..0a3393b 100644
--- a/webrtc/modules/audio_processing/audio_processing_tests.gypi
+++ b/webrtc/modules/audio_processing/audio_processing_tests.gypi
@@ -85,5 +85,22 @@
}, # click_annotate
],
}],
+ ['rtc_use_openmax_dl==1', {
+ 'targets': [
+ {
+ 'target_name': 'beamformer_test',
+ 'type': 'executable',
+ 'dependencies': [
+ '<(DEPTH)/third_party/gflags/gflags.gyp:gflags',
+ '<(webrtc_root)/modules/modules.gyp:audio_processing',
+ ],
+ 'sources': [
+ 'beamformer/beamformer_test.cc',
+ 'beamformer/pcm_utils.cc',
+ 'beamformer/pcm_utils.h',
+ ],
+ }, # beamformer_test
+ ],
+ }],
],
}
diff --git a/webrtc/modules/audio_processing/beamformer/beamformer.cc b/webrtc/modules/audio_processing/beamformer/beamformer.cc
new file mode 100644
index 0000000..4dfd3c2
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/beamformer.cc
@@ -0,0 +1,480 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#define _USE_MATH_DEFINES
+
+#include "webrtc/modules/audio_processing/beamformer/beamformer.h"
+
+#include <algorithm>
+#include <cmath>
+
+#include "webrtc/common_audio/window_generator.h"
+#include "webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h"
+
+namespace webrtc {
+namespace {
+
+// Alpha for the Kaiser Bessel Derived window.
+const float kAlpha = 1.5f;
+
+// The minimum value a postprocessing mask can take.
+const float kMaskMinimum = 0.01f;
+
+const int kFftSize = 256;
+const float kSpeedOfSoundMeterSeconds = 340;
+
+// For both target and interf angles, 0 is perpendicular to the microphone
+// array, facing forwards. The positive direction goes counterclockwise.
+// The angle at which we amplify sound.
+const float kTargetAngleRadians = 0.f;
+
+// The angle at which we suppress sound. Suppression is symmetric around 0
+// radians, so sound is suppressed at both +|kInterfAngleRadians| and
+// -|kInterfAngleRadians|. Since the beamformer is robust, this should
+// suppress sound coming from angles near +-|kInterfAngleRadians| as well.
+const float kInterfAngleRadians = static_cast<float>(M_PI) / 4.f;
+
+// When calculating the interf covariance matrix, this is the weight for
+// the weighted average between the uniform covariance matrix and the angled
+// covariance matrix.
+// Rpsi = Rpsi_angled * kBalance + Rpsi_uniform * (1 - kBalance)
+const float kBalance = 0.2f;
+
+const int kNumFreqBins = kFftSize / 2 + 1;
+
+// TODO(claguna): need comment here.
+const float kBeamwidthConstant = 0.00001f;
+
+// Width of the boxcar.
+const float kBoxcarHalfWidth = 0.001f;
+
+// We put a gap in the covariance matrix where we expect the target to come
+// from. Warning: This must be very small, ex. < 0.01, because otherwise it can
+// cause the covariance matrix not to be positive semidefinite, and we require
+// that our covariance matrices are positive semidefinite.
+const float kCovUniformGapHalfWidth = 0.001f;
+
+// How many blocks of past masks (including the current block) we save. Saved
+// masks are used for postprocessing such as removing musical noise.
+const int kNumberSavedPostfilterMasks = 2;
+
+// Lower bound on gain decay.
+const float kHalfLifeSeconds = 0.05f;
+
+// The average mask is computed from masks in this mid-frequency range.
+const int kMidFrequnecyLowerBoundHz = 250;
+const int kMidFrequencyUpperBoundHz = 400;
+
+const int kHighFrequnecyLowerBoundHz = 4000;
+const int kHighFrequencyUpperBoundHz = 7000;
+
+// Does conjugate(|norm_mat|) * |mat| * transpose(|norm_mat|). No extra space is
+// used; to accomplish this, we compute both multiplications in the same loop.
+float Norm(const ComplexMatrix<float>& mat,
+ const ComplexMatrix<float>& norm_mat) {
+ CHECK_EQ(norm_mat.num_rows(), 1);
+ CHECK_EQ(norm_mat.num_columns(), mat.num_rows());
+ CHECK_EQ(norm_mat.num_columns(), mat.num_columns());
+
+ complex<float> first_product = complex<float>(0.f, 0.f);
+ complex<float> second_product = complex<float>(0.f, 0.f);
+
+ const complex<float>* const* mat_els = mat.elements();
+ const complex<float>* const* norm_mat_els = norm_mat.elements();
+
+ for (int i = 0; i < norm_mat.num_columns(); ++i) {
+ for (int j = 0; j < norm_mat.num_columns(); ++j) {
+ complex<float> cur_norm_element = conj(norm_mat_els[0][j]);
+ complex<float> cur_mat_element = mat_els[j][i];
+ first_product += cur_norm_element * cur_mat_element;
+ }
+ second_product += first_product * norm_mat_els[0][i];
+ first_product = 0.f;
+ }
+ return second_product.real();
+}
+
+// Does conjugate(|lhs|) * |rhs| for row vectors |lhs| and |rhs|.
+complex<float> ConjugateDotProduct(const ComplexMatrix<float>& lhs,
+ const ComplexMatrix<float>& rhs) {
+ CHECK_EQ(lhs.num_rows(), 1);
+ CHECK_EQ(rhs.num_rows(), 1);
+ CHECK_EQ(lhs.num_columns(), rhs.num_columns());
+
+ const complex<float>* const* lhs_elements = lhs.elements();
+ const complex<float>* const* rhs_elements = rhs.elements();
+
+ complex<float> result = complex<float>(0.f, 0.f);
+ for (int i = 0; i < lhs.num_columns(); ++i) {
+ result += conj(lhs_elements[0][i]) * rhs_elements[0][i];
+ }
+
+ return result;
+}
+
+// Works for positive numbers only.
+int Round(float x) {
+ return std::floor(x + 0.5f);
+}
+
+} // namespace
+
+Beamformer::Beamformer(int chunk_size_ms,
+ int sample_rate_hz,
+ int num_input_channels,
+ float mic_spacing)
+ : chunk_length_(sample_rate_hz / (1000.f / chunk_size_ms)),
+ window_(new float[kFftSize]),
+ num_input_channels_(num_input_channels),
+ sample_rate_hz_(sample_rate_hz),
+ mic_spacing_(mic_spacing),
+ decay_threshold_(
+ pow(2, (kFftSize / -2.f) / (sample_rate_hz_ * kHalfLifeSeconds))),
+ mid_frequency_lower_bin_bound_(
+ Round(kMidFrequnecyLowerBoundHz * kFftSize / sample_rate_hz_)),
+ mid_frequency_upper_bin_bound_(
+ Round(kMidFrequencyUpperBoundHz * kFftSize / sample_rate_hz_)),
+ high_frequency_lower_bin_bound_(
+ Round(kHighFrequnecyLowerBoundHz * kFftSize / sample_rate_hz_)),
+ high_frequency_upper_bin_bound_(
+ Round(kHighFrequencyUpperBoundHz * kFftSize / sample_rate_hz_)),
+ current_block_ix_(0),
+ previous_block_ix_(-1),
+ postfilter_masks_(new MatrixF[kNumberSavedPostfilterMasks]),
+ delay_sum_masks_(new ComplexMatrixF[kNumFreqBins]),
+ target_cov_mats_(new ComplexMatrixF[kNumFreqBins]),
+ interf_cov_mats_(new ComplexMatrixF[kNumFreqBins]),
+ reflected_interf_cov_mats_(new ComplexMatrixF[kNumFreqBins]),
+ mask_thresholds_(new float[kNumFreqBins]),
+ wave_numbers_(new float[kNumFreqBins]),
+ rxiws_(new float[kNumFreqBins]),
+ rpsiws_(new float[kNumFreqBins]),
+ reflected_rpsiws_(new float[kNumFreqBins]) {
+ DCHECK_LE(mid_frequency_upper_bin_bound_, kNumFreqBins);
+ DCHECK_LT(mid_frequency_lower_bin_bound_, mid_frequency_upper_bin_bound_);
+ DCHECK_LE(high_frequency_upper_bin_bound_, kNumFreqBins);
+ DCHECK_LT(high_frequency_lower_bin_bound_, high_frequency_upper_bin_bound_);
+
+ WindowGenerator::KaiserBesselDerived(kAlpha, kFftSize, window_.get());
+ lapped_transform_.reset(new LappedTransform(num_input_channels_,
+ 1,
+ chunk_length_,
+ window_.get(),
+ kFftSize,
+ kFftSize / 2,
+ this));
+
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ float freq_hz = (static_cast<float>(i) / kFftSize) * sample_rate_hz_;
+ wave_numbers_[i] = 2 * M_PI * freq_hz / kSpeedOfSoundMeterSeconds;
+ }
+
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ mask_thresholds_[i] = num_input_channels_ * num_input_channels_ *
+ kBeamwidthConstant * wave_numbers_[i] *
+ wave_numbers_[i];
+ }
+
+ // Init all nonadaptive values before looping through the frames.
+ InitDelaySumMasks();
+ InitTargetCovMats();
+ InitInterfCovMats();
+
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ rxiws_[i] = Norm(target_cov_mats_[i], delay_sum_masks_[i]);
+ }
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ rpsiws_[i] = Norm(interf_cov_mats_[i], delay_sum_masks_[i]);
+ }
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ reflected_rpsiws_[i] =
+ Norm(reflected_interf_cov_mats_[i], delay_sum_masks_[i]);
+ }
+ for (int i = 0; i < kNumberSavedPostfilterMasks; ++i) {
+ postfilter_masks_[i].Resize(1, kNumFreqBins);
+ }
+}
+
+void Beamformer::InitDelaySumMasks() {
+ float sin_target = sin(kTargetAngleRadians);
+ for (int f_ix = 0; f_ix < kNumFreqBins; ++f_ix) {
+ delay_sum_masks_[f_ix].Resize(1, num_input_channels_);
+ CovarianceMatrixGenerator::PhaseAlignmentMasks(f_ix,
+ kFftSize,
+ sample_rate_hz_,
+ kSpeedOfSoundMeterSeconds,
+ mic_spacing_,
+ num_input_channels_,
+ sin_target,
+ &delay_sum_masks_[f_ix]);
+
+ complex_f norm_factor = sqrt(
+ ConjugateDotProduct(delay_sum_masks_[f_ix], delay_sum_masks_[f_ix]));
+ delay_sum_masks_[f_ix].Scale(1.f / norm_factor);
+ }
+}
+
+void Beamformer::InitTargetCovMats() {
+ target_cov_mats_[0].Resize(num_input_channels_, num_input_channels_);
+ CovarianceMatrixGenerator::DCCovarianceMatrix(
+ num_input_channels_, kBoxcarHalfWidth, &target_cov_mats_[0]);
+
+ complex_f normalization_factor = target_cov_mats_[0].Trace();
+ target_cov_mats_[0].Scale(1.f / normalization_factor);
+
+ for (int i = 1; i < kNumFreqBins; ++i) {
+ float wave_number = wave_numbers_[i];
+
+ target_cov_mats_[i].Resize(num_input_channels_, num_input_channels_);
+ CovarianceMatrixGenerator::Boxcar(wave_number,
+ num_input_channels_,
+ mic_spacing_,
+ kBoxcarHalfWidth,
+ &target_cov_mats_[i]);
+
+ complex_f normalization_factor = target_cov_mats_[i].Trace();
+ target_cov_mats_[i].Scale(1.f / normalization_factor);
+ }
+}
+
+void Beamformer::InitInterfCovMats() {
+ interf_cov_mats_[0].Resize(num_input_channels_, num_input_channels_);
+ CovarianceMatrixGenerator::DCCovarianceMatrix(
+ num_input_channels_, kCovUniformGapHalfWidth, &interf_cov_mats_[0]);
+
+ complex_f normalization_factor = interf_cov_mats_[0].Trace();
+ interf_cov_mats_[0].Scale(1.f / normalization_factor);
+
+ for (int i = 1; i < kNumFreqBins; ++i) {
+ float wave_number = wave_numbers_[i];
+
+ interf_cov_mats_[i].Resize(num_input_channels_, num_input_channels_);
+ ComplexMatrixF uniform_cov_mat(num_input_channels_, num_input_channels_);
+ ComplexMatrixF angled_cov_mat(num_input_channels_, num_input_channels_);
+
+ CovarianceMatrixGenerator::GappedUniformCovarianceMatrix(
+ wave_number,
+ num_input_channels_,
+ mic_spacing_,
+ kCovUniformGapHalfWidth,
+ &uniform_cov_mat);
+
+ CovarianceMatrixGenerator::AngledCovarianceMatrix(kSpeedOfSoundMeterSeconds,
+ kInterfAngleRadians,
+ i,
+ kFftSize,
+ kNumFreqBins,
+ sample_rate_hz_,
+ num_input_channels_,
+ mic_spacing_,
+ &angled_cov_mat);
+ // Normalize matrices before averaging them.
+ complex_f normalization_factor = uniform_cov_mat.Trace();
+ uniform_cov_mat.Scale(1.f / normalization_factor);
+ normalization_factor = angled_cov_mat.Trace();
+ angled_cov_mat.Scale(1.f / normalization_factor);
+
+ // Average matrices.
+ uniform_cov_mat.Scale(1 - kBalance);
+ angled_cov_mat.Scale(kBalance);
+ interf_cov_mats_[i].Add(uniform_cov_mat, angled_cov_mat);
+ }
+
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ reflected_interf_cov_mats_[i].PointwiseConjugate(interf_cov_mats_[i]);
+ }
+}
+
+void Beamformer::ProcessChunk(const float* const* input,
+ const float* const* high_pass_split_input,
+ int num_input_channels,
+ int num_frames_per_band,
+ float* const* output,
+ float* const* high_pass_split_output) {
+ CHECK_EQ(num_input_channels, num_input_channels_);
+ CHECK_EQ(num_frames_per_band, chunk_length_);
+
+ num_blocks_in_this_chunk_ = 0;
+ float old_high_pass_mask = high_pass_postfilter_mask_;
+ high_pass_postfilter_mask_ = 0.f;
+ high_pass_exists_ = high_pass_split_input != NULL;
+ lapped_transform_->ProcessChunk(input, output);
+
+ // Apply delay and sum and postfilter in the time domain. WARNING: only works
+ // because delay-and-sum is not frequency dependent.
+ if (high_pass_exists_) {
+ high_pass_postfilter_mask_ /= num_blocks_in_this_chunk_;
+
+ if (previous_block_ix_ == -1) {
+ old_high_pass_mask = high_pass_postfilter_mask_;
+ }
+
+ // Ramp up/down for smoothing. 1 mask per 10ms results in audible
+ // discontinuities.
+ float ramp_inc =
+ (high_pass_postfilter_mask_ - old_high_pass_mask) / num_frames_per_band;
+ for (int i = 0; i < num_frames_per_band; ++i) {
+ old_high_pass_mask += ramp_inc;
+
+ // Applying the delay and sum (at zero degrees, this is equivalent to
+ // averaging).
+ float sum = 0.f;
+ for (int j = 0; j < num_input_channels; ++j) {
+ sum += high_pass_split_input[j][i];
+ }
+ high_pass_split_output[0][i] =
+ sum / num_input_channels * old_high_pass_mask;
+ }
+ }
+}
+
+void Beamformer::ProcessAudioBlock(const complex_f* const* input,
+ int num_input_channels,
+ int num_freq_bins,
+ int num_output_channels,
+ complex_f* const* output) {
+ CHECK_EQ(num_freq_bins, kNumFreqBins);
+ CHECK_EQ(num_input_channels, num_input_channels_);
+ CHECK_EQ(num_output_channels, 1);
+
+ float* mask_data = postfilter_masks_[current_block_ix_].elements()[0];
+
+ // Calculating the postfilter masks. Note that we need two for each
+ // frequency bin to account for the positive and negative interferer
+ // angle.
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ eig_m_.CopyFromColumn(input, i, num_input_channels_);
+ float eig_m_norm_factor =
+ std::sqrt(ConjugateDotProduct(eig_m_, eig_m_)).real();
+ if (eig_m_norm_factor != 0.f) {
+ eig_m_.Scale(1.f / eig_m_norm_factor);
+ }
+
+ float rxim = Norm(target_cov_mats_[i], eig_m_);
+ float ratio_rxiw_rxim = 0.f;
+ if (rxim != 0.f) {
+ ratio_rxiw_rxim = rxiws_[i] / rxim;
+ }
+
+ complex_f rmw = abs(ConjugateDotProduct(delay_sum_masks_[i], eig_m_));
+ rmw *= rmw;
+ float rmw_r = rmw.real();
+
+ mask_data[i] = CalculatePostfilterMask(interf_cov_mats_[i],
+ rpsiws_[i],
+ ratio_rxiw_rxim,
+ rmw_r,
+ mask_thresholds_[i]);
+
+ mask_data[i] *= CalculatePostfilterMask(reflected_interf_cov_mats_[i],
+ reflected_rpsiws_[i],
+ ratio_rxiw_rxim,
+ rmw_r,
+ mask_thresholds_[i]);
+ }
+
+ // Can't access block_index - 1 on the first block.
+ if (previous_block_ix_ >= 0) {
+ ApplyDecay();
+ }
+
+ ApplyLowFrequencyCorrection();
+
+ if (high_pass_exists_) {
+ CalculateHighFrequencyMask();
+ }
+
+ ApplyMasks(input, output);
+
+ previous_block_ix_ = current_block_ix_;
+ current_block_ix_ = (current_block_ix_ + 1) % kNumberSavedPostfilterMasks;
+ num_blocks_in_this_chunk_++;
+}
+
+float Beamformer::CalculatePostfilterMask(const ComplexMatrixF& interf_cov_mat,
+ float rpsiw,
+ float ratio_rxiw_rxim,
+ float rmw_r,
+ float mask_threshold) {
+ float rpsim = Norm(interf_cov_mat, eig_m_);
+
+ // Find lambda.
+ float ratio = rpsiw / rpsim;
+ float numerator = rmw_r - ratio;
+ float denominator = ratio_rxiw_rxim - ratio;
+
+ float mask = 1.f;
+ if (denominator > mask_threshold) {
+ float lambda = numerator / denominator;
+ mask = std::max(lambda * ratio_rxiw_rxim / rmw_r, kMaskMinimum);
+ }
+ return mask;
+}
+
+void Beamformer::ApplyMasks(const complex_f* const* input,
+ complex_f* const* output) {
+ complex_f* output_channel = output[0];
+ const float* postfilter_mask_els =
+ postfilter_masks_[current_block_ix_].elements()[0];
+ for (int f_ix = 0; f_ix < kNumFreqBins; ++f_ix) {
+ output_channel[f_ix] = complex_f(0.f, 0.f);
+
+ const complex_f* delay_sum_mask_els = delay_sum_masks_[f_ix].elements()[0];
+ for (int c_ix = 0; c_ix < num_input_channels_; ++c_ix) {
+ output_channel[f_ix] += input[c_ix][f_ix] * delay_sum_mask_els[c_ix];
+ }
+
+ output_channel[f_ix] *= postfilter_mask_els[f_ix];
+ }
+}
+
+void Beamformer::ApplyDecay() {
+ float* current_mask_els = postfilter_masks_[current_block_ix_].elements()[0];
+ const float* previous_block_els =
+ postfilter_masks_[previous_block_ix_].elements()[0];
+ for (int i = 0; i < kNumFreqBins; ++i) {
+ current_mask_els[i] =
+ std::max(current_mask_els[i], previous_block_els[i] * decay_threshold_);
+ }
+}
+
+void Beamformer::ApplyLowFrequencyCorrection() {
+ float low_frequency_mask = 0.f;
+ float* mask_els = postfilter_masks_[current_block_ix_].elements()[0];
+ for (int i = mid_frequency_lower_bin_bound_;
+ i <= mid_frequency_upper_bin_bound_;
+ ++i) {
+ low_frequency_mask += mask_els[i];
+ }
+
+ low_frequency_mask /=
+ mid_frequency_upper_bin_bound_ - mid_frequency_lower_bin_bound_ + 1;
+
+ for (int i = 0; i < mid_frequency_lower_bin_bound_; ++i) {
+ mask_els[i] = low_frequency_mask;
+ }
+}
+
+void Beamformer::CalculateHighFrequencyMask() {
+ float high_pass_mask = 0.f;
+ float* mask_els = postfilter_masks_[current_block_ix_].elements()[0];
+ for (int i = high_frequency_lower_bin_bound_;
+ i <= high_frequency_upper_bin_bound_;
+ ++i) {
+ high_pass_mask += mask_els[i];
+ }
+
+ high_pass_mask /=
+ high_frequency_upper_bin_bound_ - high_frequency_lower_bin_bound_ + 1;
+
+ high_pass_postfilter_mask_ += high_pass_mask;
+}
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/beamformer.h b/webrtc/modules/audio_processing/beamformer/beamformer.h
new file mode 100644
index 0000000..31b3721
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/beamformer.h
@@ -0,0 +1,152 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_BEAMFORMER_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_BEAMFORMER_H_
+
+#include "webrtc/common_audio/lapped_transform.h"
+#include "webrtc/modules/audio_processing/beamformer/complex_matrix.h"
+
+namespace webrtc {
+
+// Enhances sound sources coming directly in front of a uniform linear array
+// and suppresses sound sources coming from all other directions. Operates on
+// multichannel signals and produces single-channel output.
+//
+// The implemented nonlinear postfilter algorithm taken from "A Robust Nonlinear
+// Beamforming Postprocessor" by Bastiaan Kleijn.
+//
+// TODO: Target angle assumed to be 0. Parameterize target angle.
+class Beamformer : public LappedTransform::Callback {
+ public:
+ Beamformer(int chunk_size_ms,
+ // Sample rate corresponds to the lower band.
+ int sample_rate_hz,
+ int num_input_channels,
+ // Microphone spacing in meters.
+ float mic_spacing);
+
+ // Process one time-domain chunk of audio. The audio can be separated into
+ // two signals by frequency, with the higher half passed in as the second
+ // parameter. Use NULL for |high_pass_split_input| if you only have one
+ // audio signal. The number of frames and channels must correspond to the
+ // ctor parameters. The same signal can be passed in as |input| and |output|.
+ void ProcessChunk(const float* const* input,
+ const float* const* high_pass_split_input,
+ int num_input_channels,
+ int num_frames_per_band,
+ float* const* output,
+ float* const* high_pass_split_output);
+
+ protected:
+ // Process one frequency-domain block of audio. This is where the fun
+ // happens. Implements LappedTransform::Callback.
+ void ProcessAudioBlock(const complex<float>* const* input,
+ int num_input_channels,
+ int num_freq_bins,
+ int num_output_channels,
+ complex<float>* const* output);
+
+ private:
+ typedef Matrix<float> MatrixF;
+ typedef ComplexMatrix<float> ComplexMatrixF;
+ typedef complex<float> complex_f;
+
+ void InitDelaySumMasks();
+ void InitTargetCovMats(); // TODO: Make this depend on target angle.
+ void InitInterfCovMats();
+
+ // An implementation of equation 18, which calculates postfilter masks that,
+ // when applied, minimize the mean-square error of our estimation of the
+ // desired signal. A sub-task is to calculate lambda, which is solved via
+ // equation 13.
+ float CalculatePostfilterMask(const ComplexMatrixF& interf_cov_mat,
+ float rpsiw,
+ float ratio_rxiw_rxim,
+ float rmxi_r,
+ float mask_threshold);
+
+ // Prevents the postfilter masks from degenerating too quickly (a cause of
+ // musical noise).
+ void ApplyDecay();
+
+ // The postfilter masks are unreliable at low frequencies. Calculates a better
+ // mask by averaging mid-low frequency values.
+ void ApplyLowFrequencyCorrection();
+
+ // Postfilter masks are also unreliable at high frequencies. Average mid-high
+ // frequency masks to calculate a single mask per block which can be applied
+ // in the time-domain. Further, we average these block-masks over a chunk,
+ // resulting in one postfilter mask per audio chunk. This allows us to skip
+ // both transforming and blocking the high-frequency signal.
+ void CalculateHighFrequencyMask();
+
+ // Applies both sets of masks to |input| and store in |output|.
+ void ApplyMasks(const complex_f* const* input, complex_f* const* output);
+
+ // Deals with the fft transform and blocking.
+ const int chunk_length_;
+ scoped_ptr<LappedTransform> lapped_transform_;
+ scoped_ptr<float[]> window_;
+
+ // Parameters exposed to the user.
+ const int num_input_channels_;
+ const int sample_rate_hz_;
+ const float mic_spacing_;
+
+ // Calculated based on user-input and constants in the .cc file.
+ const float decay_threshold_;
+ const int mid_frequency_lower_bin_bound_;
+ const int mid_frequency_upper_bin_bound_;
+ const int high_frequency_lower_bin_bound_;
+ const int high_frequency_upper_bin_bound_;
+
+ // Indices into |postfilter_masks_|.
+ int current_block_ix_;
+ int previous_block_ix_;
+
+ // Old masks are saved in this ring buffer for smoothing. Array of length
+ // |kNumberSavedMasks| matrix of size 1 x |kNumFreqBins|.
+ scoped_ptr<MatrixF[]> postfilter_masks_;
+
+ // Array of length |kNumFreqBins|, Matrix of size |1| x |num_channels_|.
+ scoped_ptr<ComplexMatrixF[]> delay_sum_masks_;
+
+ // Array of length |kNumFreqBins|, Matrix of size |num_input_channels_| x
+ // |num_input_channels_|.
+ scoped_ptr<ComplexMatrixF[]> target_cov_mats_;
+
+ // Array of length |kNumFreqBins|, Matrix of size |num_input_channels_| x
+ // |num_input_channels_|.
+ scoped_ptr<ComplexMatrixF[]> interf_cov_mats_;
+ scoped_ptr<ComplexMatrixF[]> reflected_interf_cov_mats_;
+
+ // Of length |kNumFreqBins|.
+ scoped_ptr<float[]> mask_thresholds_;
+ scoped_ptr<float[]> wave_numbers_;
+
+ // Preallocated for ProcessAudioBlock()
+ // Of length |kNumFreqBins|.
+ scoped_ptr<float[]> rxiws_;
+ scoped_ptr<float[]> rpsiws_;
+ scoped_ptr<float[]> reflected_rpsiws_;
+
+ // The microphone normalization factor.
+ ComplexMatrixF eig_m_;
+
+ // For processing the high-frequency input signal.
+ bool high_pass_exists_;
+ int num_blocks_in_this_chunk_;
+ float high_pass_postfilter_mask_;
+};
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_BEAMFORMER_H_
diff --git a/webrtc/modules/audio_processing/beamformer/beamformer_test.cc b/webrtc/modules/audio_processing/beamformer/beamformer_test.cc
new file mode 100644
index 0000000..e202506
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/beamformer_test.cc
@@ -0,0 +1,83 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#include <iostream>
+
+#include "gflags/gflags.h"
+#include "webrtc/modules/audio_processing/beamformer/beamformer.h"
+#include "webrtc/modules/audio_processing/beamformer/pcm_utils.h"
+
+DEFINE_int32(sample_rate,
+ 48000,
+ "The sample rate of the input file. The output"
+ "file will be of the same sample rate.");
+DEFINE_int32(num_input_channels,
+ 2,
+ "The number of channels in the input file.");
+DEFINE_double(mic_spacing,
+ 0.05,
+ "The spacing between microphones on the chromebook which "
+ "recorded the input file.");
+DEFINE_string(input_file_path,
+ "input.wav",
+ "The absolute path to the input file.");
+DEFINE_string(output_file_path,
+ "beamformer_test_output.wav",
+ "The absolute path to the output file.");
+
+using webrtc::ChannelBuffer;
+
+int main(int argc, char* argv[]) {
+ google::ParseCommandLineFlags(&argc, &argv, true);
+
+ const float kChunkTimeMilliseconds = 10;
+ const int kChunkSize = FLAGS_sample_rate / (1000.f / kChunkTimeMilliseconds);
+ const int kInputSamplesPerChunk = kChunkSize * FLAGS_num_input_channels;
+
+ ChannelBuffer<float> captured_audio_cb(kChunkSize, FLAGS_num_input_channels);
+
+ FILE* read_file = fopen(FLAGS_input_file_path.c_str(), "rb");
+ if (!read_file) {
+ std::cerr << "Input file '" << FLAGS_input_file_path << "' not found."
+ << std::endl;
+ return -1;
+ }
+
+ // Skipping the .wav header. TODO: Add .wav header parsing.
+ fseek(read_file, 44, SEEK_SET);
+
+ FILE* write_file = fopen(FLAGS_output_file_path.c_str(), "wb");
+ webrtc::Beamformer bf(kChunkTimeMilliseconds,
+ FLAGS_sample_rate,
+ FLAGS_num_input_channels,
+ FLAGS_mic_spacing);
+ while (true) {
+ size_t samples_read = webrtc::PcmReadToFloat(read_file,
+ kInputSamplesPerChunk,
+ FLAGS_num_input_channels,
+ captured_audio_cb.channels());
+
+ if (static_cast<int>(samples_read) != kInputSamplesPerChunk) {
+ break;
+ }
+
+ bf.ProcessChunk(captured_audio_cb.channels(),
+ NULL,
+ FLAGS_num_input_channels,
+ kChunkSize,
+ captured_audio_cb.channels(),
+ NULL);
+ webrtc::PcmWriteFromFloat(
+ write_file, kChunkSize, 1, captured_audio_cb.channels());
+ }
+ fclose(read_file);
+ fclose(write_file);
+ return 0;
+}
diff --git a/webrtc/modules/audio_processing/beamformer/complex_matrix.h b/webrtc/modules/audio_processing/beamformer/complex_matrix.h
new file mode 100644
index 0000000..0e7ee14
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/complex_matrix.h
@@ -0,0 +1,98 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_COMPLEX_MATRIX_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_COMPLEX_MATRIX_H_
+
+#include <complex>
+
+#include "webrtc/base/checks.h"
+#include "webrtc/modules/audio_processing/channel_buffer.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix.h"
+#include "webrtc/system_wrappers/interface/scoped_ptr.h"
+
+namespace webrtc {
+
+using std::complex;
+
+// An extension of Matrix for operations that only work on a complex type.
+template <typename T>
+class ComplexMatrix : public Matrix<complex<T> > {
+ public:
+ ComplexMatrix() : Matrix<complex<T> >() {}
+
+ ComplexMatrix(int num_rows, int num_columns)
+ : Matrix<complex<T> >(num_rows, num_columns) {}
+
+ ComplexMatrix(const complex<T>* data, int num_rows, int num_columns)
+ : Matrix<complex<T> >(data, num_rows, num_columns) {}
+
+ // Complex Matrix operations.
+ ComplexMatrix& PointwiseConjugate() {
+ complex<T>* const data = this->data();
+ size_t size = this->num_rows() * this->num_columns();
+ for (size_t i = 0; i < size; ++i) {
+ data[i] = conj(data[i]);
+ }
+
+ return *this;
+ }
+
+ ComplexMatrix& PointwiseConjugate(const ComplexMatrix& operand) {
+ this->CopyFrom(operand);
+ return PointwiseConjugate();
+ }
+
+ ComplexMatrix& ConjugateTranspose() {
+ this->CopyDataToScratch();
+ int num_rows = this->num_rows();
+ this->SetNumRows(this->num_columns());
+ this->SetNumColumns(num_rows);
+ this->Resize();
+ return ConjugateTranspose(this->scratch_elements());
+ }
+
+ ComplexMatrix& ConjugateTranspose(const ComplexMatrix& operand) {
+ CHECK_EQ(operand.num_rows(), this->num_columns());
+ CHECK_EQ(operand.num_columns(), this->num_rows());
+ return ConjugateTranspose(operand.elements());
+ }
+
+ ComplexMatrix& ZeroImag() {
+ complex<T>* const data = this->data();
+ size_t size = this->num_rows() * this->num_columns();
+ for (size_t i = 0; i < size; ++i) {
+ data[i] = complex<T>(data[i].real(), 0);
+ }
+
+ return *this;
+ }
+
+ ComplexMatrix& ZeroImag(const ComplexMatrix& operand) {
+ this->CopyFrom(operand);
+ return ZeroImag();
+ }
+
+ private:
+ ComplexMatrix& ConjugateTranspose(const complex<T>* const* src) {
+ complex<T>* const* elements = this->elements();
+ for (int i = 0; i < this->num_rows(); ++i) {
+ for (int j = 0; j < this->num_columns(); ++j) {
+ elements[i][j] = conj(src[j][i]);
+ }
+ }
+
+ return *this;
+ }
+};
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_COMPLEX_MATRIX_H_
diff --git a/webrtc/modules/audio_processing/beamformer/complex_matrix_unittest.cc b/webrtc/modules/audio_processing/beamformer/complex_matrix_unittest.cc
new file mode 100644
index 0000000..63c8c4f
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/complex_matrix_unittest.cc
@@ -0,0 +1,102 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#include "testing/gtest/include/gtest/gtest.h"
+#include "webrtc/modules/audio_processing/beamformer/complex_matrix.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h"
+
+namespace webrtc {
+
+TEST(ComplexMatrixTest, TestPointwiseConjugate) {
+ const int kNumRows = 2;
+ const int kNumCols = 4;
+
+ const complex<float> kValuesInitial[kNumRows][kNumCols] = {
+ {complex<float>(1.1f, 1.1f), complex<float>(2.2f, -2.2f),
+ complex<float>(3.3f, 3.3f), complex<float>(4.4f, -4.4f)},
+ {complex<float>(5.5f, 5.5f), complex<float>(6.6f, -6.6f),
+ complex<float>(7.7f, 7.7f), complex<float>(8.8f, -8.8f)}};
+
+ const complex<float> kValuesExpected[kNumRows][kNumCols] = {
+ {complex<float>(1.1f, -1.1f), complex<float>(2.2f, 2.2f),
+ complex<float>(3.3f, -3.3f), complex<float>(4.4f, 4.4f)},
+ {complex<float>(5.5f, -5.5f), complex<float>(6.6f, 6.6f),
+ complex<float>(7.7f, -7.7f), complex<float>(8.8f, 8.8f)}};
+
+ ComplexMatrix<float> initial_mat(*kValuesInitial, kNumRows, kNumCols);
+ ComplexMatrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ ComplexMatrix<float> actual_result(kNumRows, kNumCols);
+
+ actual_result.PointwiseConjugate(initial_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result,
+ actual_result);
+
+ initial_mat.PointwiseConjugate();
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(initial_mat,
+ actual_result);
+}
+
+TEST(ComplexMatrixTest, TestConjugateTranspose) {
+ const int kNumInitialRows = 2;
+ const int kNumInitialCols = 4;
+ const int kNumResultRows = 4;
+ const int kNumResultCols = 2;
+
+ const complex<float> kValuesInitial[kNumInitialRows][kNumInitialCols] = {
+ {complex<float>(1.1f, 1.1f), complex<float>(2.2f, 2.2f),
+ complex<float>(3.3f, 3.3f), complex<float>(4.4f, 4.4f)},
+ {complex<float>(5.5f, 5.5f), complex<float>(6.6f, 6.6f),
+ complex<float>(7.7f, 7.7f), complex<float>(8.8f, 8.8f)}};
+
+ const complex<float> kValuesExpected[kNumResultRows][kNumResultCols] = {
+ {complex<float>(1.1f, -1.1f), complex<float>(5.5f, -5.5f)},
+ {complex<float>(2.2f, -2.2f), complex<float>(6.6f, -6.6f)},
+ {complex<float>(3.3f, -3.3f), complex<float>(7.7f, -7.7f)},
+ {complex<float>(4.4f, -4.4f), complex<float>(8.8f, -8.8f)}};
+
+ ComplexMatrix<float> initial_mat(
+ *kValuesInitial, kNumInitialRows, kNumInitialCols);
+ ComplexMatrix<float> expected_result(
+ *kValuesExpected, kNumResultRows, kNumResultCols);
+ ComplexMatrix<float> actual_result(kNumResultRows, kNumResultCols);
+
+ actual_result.ConjugateTranspose(initial_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result,
+ actual_result);
+
+ initial_mat.ConjugateTranspose();
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(initial_mat,
+ actual_result);
+}
+
+TEST(ComplexMatrixTest, TestZeroImag) {
+ const int kNumRows = 2;
+ const int kNumCols = 2;
+ const complex<float> kValuesInitial[kNumRows][kNumCols] = {
+ {complex<float>(1.1f, 1.1f), complex<float>(2.2f, 2.2f)},
+ {complex<float>(3.3f, 3.3f), complex<float>(4.4f, 4.4f)}};
+ const complex<float> kValuesExpected[kNumRows][kNumCols] = {
+ {complex<float>(1.1f, 0.f), complex<float>(2.2f, 0.f)},
+ {complex<float>(3.3f, 0.f), complex<float>(4.4f, 0.f)}};
+
+ ComplexMatrix<float> initial_mat(*kValuesInitial, kNumRows, kNumCols);
+ ComplexMatrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ ComplexMatrix<float> actual_result;
+
+ actual_result.ZeroImag(initial_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result,
+ actual_result);
+
+ initial_mat.ZeroImag();
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(initial_mat,
+ actual_result);
+}
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.cc b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.cc
new file mode 100644
index 0000000..856428e
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.cc
@@ -0,0 +1,154 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#define _USE_MATH_DEFINES
+
+#include "webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h"
+
+#include <cmath>
+
+namespace {
+
+float BesselJ0(float x) {
+#if WEBRTC_WIN
+ return _j0(x);
+#else
+ return j0(x);
+#endif
+}
+
+} // namespace
+
+namespace webrtc {
+
+// Calculates the boxcar-angular desired source distribution at a given
+// wavenumber, and stores it in |mat|.
+void CovarianceMatrixGenerator::Boxcar(float wave_number,
+ int num_input_channels,
+ float mic_spacing,
+ float half_width,
+ ComplexMatrix<float>* mat) {
+ CHECK_EQ(num_input_channels, mat->num_rows());
+ CHECK_EQ(num_input_channels, mat->num_columns());
+
+ complex<float>* const* boxcar_elements = mat->elements();
+
+ for (int i = 0; i < num_input_channels; ++i) {
+ for (int j = 0; j < num_input_channels; ++j) {
+ if (i == j) {
+ boxcar_elements[i][j] = complex<float>(2.f * half_width, 0.f);
+ } else {
+ float factor = (j - i) * wave_number * mic_spacing;
+ float boxcar_real = 2.f * sin(factor * half_width) / factor;
+ boxcar_elements[i][j] = complex<float>(boxcar_real, 0.f);
+ }
+ }
+ }
+}
+
+void CovarianceMatrixGenerator::GappedUniformCovarianceMatrix(
+ float wave_number,
+ float num_input_channels,
+ float mic_spacing,
+ float gap_half_width,
+ ComplexMatrix<float>* mat) {
+ CHECK_EQ(num_input_channels, mat->num_rows());
+ CHECK_EQ(num_input_channels, mat->num_columns());
+
+ complex<float>* const* mat_els = mat->elements();
+ for (int i = 0; i < num_input_channels; ++i) {
+ for (int j = 0; j < num_input_channels; ++j) {
+ float x = (j - i) * wave_number * mic_spacing;
+ mat_els[i][j] = BesselJ0(x);
+ }
+ }
+
+ ComplexMatrix<float> boxcar_mat(num_input_channels, num_input_channels);
+ CovarianceMatrixGenerator::Boxcar(wave_number,
+ num_input_channels,
+ mic_spacing,
+ gap_half_width,
+ &boxcar_mat);
+ mat->Subtract(boxcar_mat);
+}
+
+void CovarianceMatrixGenerator::AngledCovarianceMatrix(
+ float sound_speed,
+ float angle,
+ int frequency_bin,
+ int fft_size,
+ int num_freq_bins,
+ int sample_rate,
+ int num_input_channels,
+ float mic_spacing,
+ ComplexMatrix<float>* mat) {
+ CHECK_EQ(num_input_channels, mat->num_rows());
+ CHECK_EQ(num_input_channels, mat->num_columns());
+
+ ComplexMatrix<float> interf_cov_vector(1, num_input_channels);
+ ComplexMatrix<float> interf_cov_vector_transposed(num_input_channels, 1);
+ PhaseAlignmentMasks(frequency_bin,
+ fft_size,
+ sample_rate,
+ sound_speed,
+ mic_spacing,
+ num_input_channels,
+ sin(angle),
+ &interf_cov_vector);
+ interf_cov_vector_transposed.Transpose(interf_cov_vector);
+ interf_cov_vector.PointwiseConjugate();
+ mat->Multiply(interf_cov_vector_transposed, interf_cov_vector);
+}
+
+void CovarianceMatrixGenerator::DCCovarianceMatrix(int num_input_channels,
+ float half_width,
+ ComplexMatrix<float>* mat) {
+ CHECK_EQ(num_input_channels, mat->num_rows());
+ CHECK_EQ(num_input_channels, mat->num_columns());
+
+ complex<float>* const* elements = mat->elements();
+
+ float diagonal_value = 1 - (2 * half_width);
+ for (int i = 0; i < num_input_channels; ++i) {
+ for (int j = 0; j < num_input_channels; ++j) {
+ if (i == j) {
+ elements[i][j] = complex<float>(diagonal_value, 0.f);
+ } else {
+ elements[i][j] = complex<float>(0.f, 0.f);
+ }
+ }
+ }
+}
+
+void CovarianceMatrixGenerator::PhaseAlignmentMasks(int frequency_bin,
+ int fft_size,
+ int sample_rate,
+ float sound_speed,
+ float mic_spacing,
+ int num_input_channels,
+ float sin_angle,
+ ComplexMatrix<float>* mat) {
+ CHECK_EQ(1, mat->num_rows());
+ CHECK_EQ(num_input_channels, mat->num_columns());
+
+ float freq_in_hertz =
+ (static_cast<float>(frequency_bin) / fft_size) * sample_rate;
+
+ complex<float>* const* mat_els = mat->elements();
+ for (int c_ix = 0; c_ix < num_input_channels; ++c_ix) {
+ float distance = mic_spacing * c_ix * sin_angle * -1.f;
+ float phase_shift = 2 * M_PI * distance * freq_in_hertz / sound_speed;
+
+ // Euler's formula for mat[0][c_ix] = e^(j * phase_shift).
+ mat_els[0][c_ix] = complex<float>(cos(phase_shift), sin(phase_shift));
+ }
+}
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h
new file mode 100644
index 0000000..a3bceef
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h
@@ -0,0 +1,72 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_COVARIANCE_MATRIX_GENERATOR_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_COVARIANCE_MATRIX_GENERATOR_H_
+
+#include "webrtc/modules/audio_processing/beamformer/complex_matrix.h"
+
+namespace webrtc {
+
+// Helper class for Beamformer in charge of generating covariance matrices. For
+// each function, the passed-in ComplexMatrix is expected to be of size
+// |num_input_channels| x |num_input_channels|.
+class CovarianceMatrixGenerator {
+ public:
+ // Generates the covariance matrix of the target. The boxcar implementation
+ // suppresses some high-frequency distortion caused by narrow high-frequency
+ // suppression bands turning on/off too quickly. WARNING: The target angle is
+ // assumed to be 0.
+ static void Boxcar(float wave_number,
+ int num_input_channels,
+ float mic_spacing,
+ float half_width,
+ ComplexMatrix<float>* mat);
+
+ // A uniform covariance matrix with a gap at the target location. WARNING:
+ // The target angle is assumed to be 0.
+ static void GappedUniformCovarianceMatrix(float wave_number,
+ float num_input_channels,
+ float mic_spacing,
+ float gap_half_width,
+ ComplexMatrix<float>* mat);
+
+ // The covariance matrix of a source at the given angle.
+ static void AngledCovarianceMatrix(float sound_speed,
+ float angle,
+ int frequency_bin,
+ int fft_size,
+ int num_freq_bins,
+ int sample_rate,
+ int num_input_channels,
+ float mic_spacing,
+ ComplexMatrix<float>* mat);
+
+ // A base-case covariance matrix for when the frequency is 0 Hertz.
+ static void DCCovarianceMatrix(int num_input_channels,
+ float half_width,
+ ComplexMatrix<float>* mat);
+
+ // Calculates phase shifts that, when applied to a multichannel signal and
+ // added together, cause constructive interferernce for sources located at
+ // the given angle.
+ static void PhaseAlignmentMasks(int frequency_bin,
+ int fft_size,
+ int sample_rate,
+ float sound_speed,
+ float mic_spacing,
+ int num_input_channels,
+ float sin_angle,
+ ComplexMatrix<float>* mat);
+};
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_BF_HELPERS_H_
diff --git a/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator_unittest.cc b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator_unittest.cc
new file mode 100644
index 0000000..0c3a1d5
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/covariance_matrix_generator_unittest.cc
@@ -0,0 +1,263 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#define _USE_MATH_DEFINES
+
+#include "webrtc/modules/audio_processing/beamformer/covariance_matrix_generator.h"
+
+#include <cmath>
+
+#include "testing/gtest/include/gtest/gtest.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h"
+
+namespace webrtc {
+
+using std::complex;
+
+TEST(CovarianceMatrixGeneratorTest, TestBoxcar) {
+ const float kWaveNumber = 10.3861f;
+ const int kNumberMics = 3;
+ const float kMicSpacing = 0.23f;
+ const float kHalfWidth = 0.001f;
+ const float kTolerance = 0.0001f;
+
+ ComplexMatrix<float> actual_boxcar(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::Boxcar(
+ kWaveNumber, kNumberMics, kMicSpacing, kHalfWidth, &actual_boxcar);
+
+ complex<float>* const* actual_els = actual_boxcar.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].real(), 0.002f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].real(), 0.002f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].imag(), 0.f, kTolerance);
+}
+
+// Test Boxcar with large numbers to generate matrix with varying
+// numbers.
+TEST(CovarianceMatrixGeneratorTest, TestBoxcarLargeNumbers) {
+ const float kWaveNumber = 10.3861f;
+ const int kNumberMics = 3;
+ const float kMicSpacing = 3.f;
+ const float kHalfWidth = 0.1f;
+ const float kTolerance = 0.0001f;
+
+ ComplexMatrix<float> actual_boxcar(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::Boxcar(
+ kWaveNumber, kNumberMics, kMicSpacing, kHalfWidth, &actual_boxcar);
+
+ complex<float>* const* actual_els = actual_boxcar.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 0.2f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].real(), -0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 0.2f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].real(), 0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].real(), -0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].real(), 0.0017f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].real(), 0.2f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].imag(), 0.f, kTolerance);
+}
+
+TEST(CovarianceMatrixGeneratorTest, TestGappedUniformCovarianceMatrix2Mics) {
+ const float kWaveNumber = 0.5775f;
+ const int kNumberMics = 2;
+ const float kMicSpacing = 0.05f;
+ const float kHalfWidth = 0.001f;
+ const float kTolerance = 0.0001f;
+
+ ComplexMatrix<float> actual_covariance_matrix(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::GappedUniformCovarianceMatrix(
+ kWaveNumber,
+ kNumberMics,
+ kMicSpacing,
+ kHalfWidth,
+ &actual_covariance_matrix);
+
+ complex<float>* const* actual_els = actual_covariance_matrix.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 0.998f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.9978f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.9978f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 0.998f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+}
+
+TEST(CovarianceMatrixGeneratorTest, TestGappedUniformCovarianceMatrix3Mics) {
+ const float kWaveNumber = 10.3861f;
+ const int kNumberMics = 3;
+ const float kMicSpacing = 0.04f;
+ const float kHalfWidth = 0.001f;
+ const float kTolerance = 0.0001f;
+
+ ComplexMatrix<float> actual_covariance_matrix(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::GappedUniformCovarianceMatrix(
+ kWaveNumber,
+ kNumberMics,
+ kMicSpacing,
+ kHalfWidth,
+ &actual_covariance_matrix);
+
+ complex<float>* const* actual_els = actual_covariance_matrix.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 0.998f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.9553f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].real(), 0.8327f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.9553f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 0.998f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].real(), 0.9553f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].real(), 0.8327f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].real(), 0.9553f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].real(), 0.998f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].imag(), 0.f, kTolerance);
+}
+
+TEST(CovarianceMatrixGeneratorTest, TestAngledCovarianceMatrix2Mics) {
+ const float kSpeedOfSound = 340;
+ const float kAngle = static_cast<float>(M_PI) / 4.f;
+ const float kFrequencyBin = 6;
+ const float kFftSize = 512;
+ const int kNumberFrequencyBins = 257;
+ const int kSampleRate = 16000;
+ const int kNumberMics = 2;
+ const float kMicSpacing = 0.04f;
+ const float kTolerance = 0.0001f;
+ ComplexMatrix<float> actual_covariance_matrix(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::AngledCovarianceMatrix(kSpeedOfSound,
+ kAngle,
+ kFrequencyBin,
+ kFftSize,
+ kNumberFrequencyBins,
+ kSampleRate,
+ kNumberMics,
+ kMicSpacing,
+ &actual_covariance_matrix);
+
+ complex<float>* const* actual_els = actual_covariance_matrix.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 1.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.9952f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.9952f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 1.f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.0978f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), -0.0978f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+}
+
+TEST(CovarianceMatrixGeneratorTest, TestAngledCovarianceMatrix3Mics) {
+ const float kSpeedOfSound = 340;
+ const float kAngle = static_cast<float>(M_PI) / 4.f;
+ const float kFrequencyBin = 9;
+ const float kFftSize = 512;
+ const int kNumberFrequencyBins = 257;
+ const int kSampleRate = 42000;
+ const int kNumberMics = 3;
+ const float kMicSpacing = 0.05f;
+ const float kTolerance = 0.0001f;
+
+ ComplexMatrix<float> actual_covariance_matrix(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::AngledCovarianceMatrix(kSpeedOfSound,
+ kAngle,
+ kFrequencyBin,
+ kFftSize,
+ kNumberFrequencyBins,
+ kSampleRate,
+ kNumberMics,
+ kMicSpacing,
+ &actual_covariance_matrix);
+
+ complex<float>* const* actual_els = actual_covariance_matrix.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 1.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.8859f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].real(), 0.5696f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.8859f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 1.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].real(), 0.8859f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].real(), 0.5696f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].real(), 0.8859f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].real(), 1.f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.4639f, kTolerance);
+ EXPECT_NEAR(actual_els[0][2].imag(), 0.8219f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), -0.4639f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][2].imag(), 0.4639f, kTolerance);
+ EXPECT_NEAR(actual_els[2][0].imag(), -0.8219f, kTolerance);
+ EXPECT_NEAR(actual_els[2][1].imag(), -0.4639f, kTolerance);
+ EXPECT_NEAR(actual_els[2][2].imag(), 0.f, kTolerance);
+}
+
+TEST(CovarianceMatrixGeneratorTest, TestDCCovarianceMatrix) {
+ const int kNumberMics = 2;
+ const float kHalfWidth = 0.01f;
+
+ ComplexMatrix<float> actual_covariance_matrix(kNumberMics, kNumberMics);
+ CovarianceMatrixGenerator::DCCovarianceMatrix(
+ kNumberMics, kHalfWidth, &actual_covariance_matrix);
+
+ complex<float>* const* actual_els = actual_covariance_matrix.elements();
+
+ EXPECT_NEAR(actual_els[0][0].real(), 0.98f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].real(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].real(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].real(), 0.98f, kTolerance);
+
+ EXPECT_NEAR(actual_els[0][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[0][1].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][0].imag(), 0.f, kTolerance);
+ EXPECT_NEAR(actual_els[1][1].imag(), 0.f, kTolerance);
+}
+
+// PhaseAlignmentMasks is tested by AngledCovarianceMatrix and by
+// InitBeamformerWeights in BeamformerUnittest.
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/matrix.h b/webrtc/modules/audio_processing/beamformer/matrix.h
new file mode 100644
index 0000000..1cb475c
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/matrix.h
@@ -0,0 +1,373 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_H_
+
+#include <algorithm>
+#include <string>
+#include <vector>
+
+#include "webrtc/base/checks.h"
+#include "webrtc/base/constructormagic.h"
+#include "webrtc/modules/audio_processing/channel_buffer.h"
+#include "webrtc/system_wrappers/interface/scoped_ptr.h"
+
+namespace {
+
+// Wrappers to get around the compiler warning resulting from the fact that
+// there's no std::sqrt overload for ints. We cast all non-complex types to
+// a double for the sqrt method.
+template <typename T>
+T sqrt_wrapper(T x) {
+ return sqrt(static_cast<double>(x));
+}
+
+template <typename S>
+std::complex<S> sqrt_wrapper(std::complex<S> x) {
+ return sqrt(x);
+}
+} // namespace
+
+namespace webrtc {
+
+// Matrix is a class for doing standard matrix operations on 2 dimensional
+// matrices of any size. Results of matrix operations are stored in the
+// calling object. Function overloads exist for both in-place (the calling
+// object is used as both an operand and the result) and out-of-place (all
+// operands are passed in as parameters) operations. If operand dimensions
+// mismatch, the program crashes. Out-of-place operations change the size of
+// the calling object, if necessary, before operating.
+//
+// 'In-place' operations that inherently change the size of the matrix (eg.
+// Transpose, Multiply on different-sized matrices) must make temporary copies
+// (|scratch_elements_| and |scratch_data_|) of existing data to complete the
+// operations.
+//
+// The data is stored contiguously. Data can be accessed internally as a flat
+// array, |data_|, or as an array of row pointers, |elements_|, but is
+// available to users only as an array of row pointers through |elements()|.
+// Memory for storage is allocated when a matrix is resized only if the new
+// size overflows capacity. Memory needed temporarily for any operations is
+// similarly resized only if the new size overflows capacity.
+//
+// If you pass in storage through the ctor, that storage is copied into the
+// matrix. TODO(claguna): albeit tricky, allow for data to be referenced
+// instead of copied, and owned by the user.
+template <typename T>
+class Matrix {
+ public:
+ Matrix() : num_rows_(0), num_columns_(0) {}
+
+ // Allocates space for the elements and initializes all values to zero.
+ Matrix(int num_rows, int num_columns)
+ : num_rows_(num_rows), num_columns_(num_columns) {
+ Resize();
+ scratch_data_.resize(num_rows_ * num_columns_);
+ scratch_elements_.resize(num_rows_);
+ }
+
+ // Copies |data| into the new Matrix.
+ Matrix(const T* data, int num_rows, int num_columns)
+ : num_rows_(num_rows), num_columns_(num_columns) {
+ CopyFrom(data, num_rows, num_columns);
+ scratch_data_.resize(num_rows_ * num_columns_);
+ scratch_elements_.resize(num_rows_);
+ }
+
+ virtual ~Matrix() {}
+
+ // Deep copy an existing matrix.
+ void CopyFrom(const Matrix& other) {
+ CopyFrom(&other.data_[0], other.num_rows_, other.num_columns_);
+ }
+
+ // Copy |data| into the Matrix. The current data is lost.
+ void CopyFrom(const T* const data, int num_rows, int num_columns) {
+ num_rows_ = num_rows;
+ num_columns_ = num_columns;
+ size_t size = num_rows_ * num_columns_;
+
+ data_.assign(data, data + size);
+ elements_.resize(num_rows_);
+ for (int i = 0; i < num_rows_; ++i) {
+ elements_[i] = &data_[i * num_columns_];
+ }
+ }
+
+ Matrix& CopyFromColumn(const T* const* src, int column_index, int num_rows) {
+ Resize(1, num_rows);
+ for (int i = 0; i < num_columns_; ++i) {
+ data_[i] = src[i][column_index];
+ }
+
+ return *this;
+ }
+
+ void Resize(int num_rows, int num_columns) {
+ num_rows_ = num_rows;
+ num_columns_ = num_columns;
+ Resize();
+ }
+
+ // Accessors and mutators.
+ int num_rows() const { return num_rows_; }
+ int num_columns() const { return num_columns_; }
+ T* const* elements() { return &elements_[0]; }
+ const T* const* elements() const { return &elements_[0]; }
+
+ T Trace() {
+ CHECK_EQ(num_rows_, num_columns_);
+
+ T trace = 0;
+ for (int i = 0; i < num_rows_; ++i) {
+ trace += elements_[i][i];
+ }
+ return trace;
+ }
+
+ // Matrix Operations. Returns *this to support method chaining.
+ Matrix& Transpose() {
+ CopyDataToScratch();
+ std::swap(num_rows_, num_columns_);
+ Resize();
+ return Transpose(scratch_elements());
+ }
+
+ Matrix& Transpose(const Matrix& operand) {
+ CHECK_EQ(operand.num_rows_, num_columns_);
+ CHECK_EQ(operand.num_columns_, num_rows_);
+
+ return Transpose(operand.elements());
+ }
+
+ template <typename S>
+ Matrix& Scale(const S& scalar) {
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] *= scalar;
+ }
+
+ return *this;
+ }
+
+ template <typename S>
+ Matrix& Scale(const Matrix& operand, const S& scalar) {
+ CopyFrom(operand);
+ return Scale(scalar);
+ }
+
+ Matrix& Add(const Matrix& operand) {
+ CHECK_EQ(num_rows_, operand.num_rows_);
+ CHECK_EQ(num_columns_, operand.num_columns_);
+
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] += operand.data_[i];
+ }
+
+ return *this;
+ }
+
+ Matrix& Add(const Matrix& lhs, const Matrix& rhs) {
+ CopyFrom(lhs);
+ return Add(rhs);
+ }
+
+ Matrix& Subtract(const Matrix& operand) {
+ CHECK_EQ(num_rows_, operand.num_rows_);
+ CHECK_EQ(num_columns_, operand.num_columns_);
+
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] -= operand.data_[i];
+ }
+
+ return *this;
+ }
+
+ Matrix& Subtract(const Matrix& lhs, const Matrix& rhs) {
+ CopyFrom(lhs);
+ return Subtract(rhs);
+ }
+
+ Matrix& PointwiseMultiply(const Matrix& operand) {
+ CHECK_EQ(num_rows_, operand.num_rows_);
+ CHECK_EQ(num_columns_, operand.num_columns_);
+
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] *= operand.data_[i];
+ }
+
+ return *this;
+ }
+
+ Matrix& PointwiseMultiply(const Matrix& lhs, const Matrix& rhs) {
+ CopyFrom(lhs);
+ return PointwiseMultiply(rhs);
+ }
+
+ Matrix& PointwiseDivide(const Matrix& operand) {
+ CHECK_EQ(num_rows_, operand.num_rows_);
+ CHECK_EQ(num_columns_, operand.num_columns_);
+
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] /= operand.data_[i];
+ }
+
+ return *this;
+ }
+
+ Matrix& PointwiseDivide(const Matrix& lhs, const Matrix& rhs) {
+ CopyFrom(lhs);
+ return PointwiseDivide(rhs);
+ }
+
+ Matrix& PointwiseSquareRoot() {
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] = sqrt_wrapper(data_[i]);
+ }
+
+ return *this;
+ }
+
+ Matrix& PointwiseSquareRoot(const Matrix& operand) {
+ CopyFrom(operand);
+ return PointwiseSquareRoot();
+ }
+
+ Matrix& PointwiseAbsoluteValue() {
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] = abs(data_[i]);
+ }
+
+ return *this;
+ }
+
+ Matrix& PointwiseAbsoluteValue(const Matrix& operand) {
+ CopyFrom(operand);
+ return PointwiseAbsoluteValue();
+ }
+
+ Matrix& PointwiseSquare() {
+ for (size_t i = 0; i < data_.size(); ++i) {
+ data_[i] *= data_[i];
+ }
+
+ return *this;
+ }
+
+ Matrix& PointwiseSquare(const Matrix& operand) {
+ CopyFrom(operand);
+ return PointwiseSquare();
+ }
+
+ Matrix& Multiply(const Matrix& lhs, const Matrix& rhs) {
+ CHECK_EQ(lhs.num_columns_, rhs.num_rows_);
+ CHECK_EQ(num_rows_, lhs.num_rows_);
+ CHECK_EQ(num_columns_, rhs.num_columns_);
+
+ return Multiply(lhs.elements(), rhs.num_rows_, rhs.elements());
+ }
+
+ Matrix& Multiply(const Matrix& rhs) {
+ CHECK_EQ(num_columns_, rhs.num_rows_);
+
+ CopyDataToScratch();
+ num_columns_ = rhs.num_columns_;
+ Resize();
+ return Multiply(scratch_elements(), rhs.num_rows_, rhs.elements());
+ }
+
+ std::string ToString() const {
+ std::ostringstream ss;
+ ss << std::endl << "Matrix" << std::endl;
+
+ for (int i = 0; i < num_rows_; ++i) {
+ for (int j = 0; j < num_columns_; ++j) {
+ ss << elements_[i][j] << " ";
+ }
+ ss << std::endl;
+ }
+ ss << std::endl;
+
+ return ss.str();
+ }
+
+ protected:
+ void SetNumRows(const int num_rows) { num_rows_ = num_rows; }
+ void SetNumColumns(const int num_columns) { num_columns_ = num_columns; }
+ T* data() { return &data_[0]; }
+ const T* data() const { return &data_[0]; }
+ const T* const* scratch_elements() const { return &scratch_elements_[0]; }
+
+ // Resize the matrix. If an increase in capacity is required, the current
+ // data is lost.
+ void Resize() {
+ size_t size = num_rows_ * num_columns_;
+ data_.resize(size);
+ elements_.resize(num_rows_);
+
+ for (int i = 0; i < num_rows_; ++i) {
+ elements_[i] = &data_[i * num_columns_];
+ }
+ }
+
+ // Copies data_ into scratch_data_ and updates scratch_elements_ accordingly.
+ void CopyDataToScratch() {
+ scratch_data_ = data_;
+ scratch_elements_.resize(num_rows_);
+
+ for (int i = 0; i < num_rows_; ++i) {
+ scratch_elements_[i] = &scratch_data_[i * num_columns_];
+ }
+ }
+
+ private:
+ int num_rows_;
+ int num_columns_;
+ std::vector<T> data_;
+ std::vector<T*> elements_;
+
+ // Stores temporary copies of |data_| and |elements_| for in-place operations
+ // where referring to original data is necessary.
+ std::vector<T> scratch_data_;
+ std::vector<T*> scratch_elements_;
+
+ // Helpers for Transpose and Multiply operations that unify in-place and
+ // out-of-place solutions.
+ Matrix& Transpose(const T* const* src) {
+ for (int i = 0; i < num_rows_; ++i) {
+ for (int j = 0; j < num_columns_; ++j) {
+ elements_[i][j] = src[j][i];
+ }
+ }
+
+ return *this;
+ }
+
+ Matrix& Multiply(const T* const* lhs, int num_rows_rhs, const T* const* rhs) {
+ for (int row = 0; row < num_rows_; ++row) {
+ for (int col = 0; col < num_columns_; ++col) {
+ T cur_element = 0;
+ for (int i = 0; i < num_rows_rhs; ++i) {
+ cur_element += lhs[row][i] * rhs[i][col];
+ }
+
+ elements_[row][col] = cur_element;
+ }
+ }
+
+ return *this;
+ }
+
+ DISALLOW_COPY_AND_ASSIGN(Matrix);
+};
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_H_
diff --git a/webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h b/webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h
new file mode 100644
index 0000000..7c58670
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h
@@ -0,0 +1,102 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_TEST_HELPERS_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_TEST_HELPERS_H_
+
+#include "testing/gtest/include/gtest/gtest.h"
+#include "webrtc/modules/audio_processing/beamformer/complex_matrix.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix.h"
+
+namespace {
+const float kTolerance = 0.001f;
+}
+
+namespace webrtc {
+
+using std::complex;
+
+// Functions used in both matrix_unittest and complex_matrix_unittest.
+class MatrixTestHelpers {
+ public:
+ template <typename T>
+ static void ValidateMatrixEquality(const Matrix<T>& expected,
+ const Matrix<T>& actual) {
+ EXPECT_EQ(expected.num_rows(), actual.num_rows());
+ EXPECT_EQ(expected.num_columns(), actual.num_columns());
+
+ const T* const* expected_elements = expected.elements();
+ const T* const* actual_elements = actual.elements();
+ for (int i = 0; i < expected.num_rows(); ++i) {
+ for (int j = 0; j < expected.num_columns(); ++j) {
+ EXPECT_EQ(expected_elements[i][j], actual_elements[i][j]);
+ }
+ }
+ }
+
+ static void ValidateMatrixEqualityFloat(const Matrix<float>& expected,
+ const Matrix<float>& actual) {
+ EXPECT_EQ(expected.num_rows(), actual.num_rows());
+ EXPECT_EQ(expected.num_columns(), actual.num_columns());
+
+ const float* const* expected_elements = expected.elements();
+ const float* const* actual_elements = actual.elements();
+ for (int i = 0; i < expected.num_rows(); ++i) {
+ for (int j = 0; j < expected.num_columns(); ++j) {
+ EXPECT_NEAR(expected_elements[i][j], actual_elements[i][j], kTolerance);
+ }
+ }
+ }
+
+ static void ValidateMatrixEqualityComplexFloat(
+ const Matrix<complex<float> >& expected,
+ const Matrix<complex<float> >& actual) {
+ EXPECT_EQ(expected.num_rows(), actual.num_rows());
+ EXPECT_EQ(expected.num_columns(), actual.num_columns());
+
+ const complex<float>* const* expected_elements = expected.elements();
+ const complex<float>* const* actual_elements = actual.elements();
+ for (int i = 0; i < expected.num_rows(); ++i) {
+ for (int j = 0; j < expected.num_columns(); ++j) {
+ EXPECT_NEAR(expected_elements[i][j].real(),
+ actual_elements[i][j].real(),
+ kTolerance);
+ EXPECT_NEAR(expected_elements[i][j].imag(),
+ actual_elements[i][j].imag(),
+ kTolerance);
+ }
+ }
+ }
+
+ static void ValidateMatrixNearEqualityComplexFloat(
+ const Matrix<complex<float> >& expected,
+ const Matrix<complex<float> >& actual,
+ float tolerance) {
+ EXPECT_EQ(expected.num_rows(), actual.num_rows());
+ EXPECT_EQ(expected.num_columns(), actual.num_columns());
+
+ const complex<float>* const* expected_elements = expected.elements();
+ const complex<float>* const* actual_elements = actual.elements();
+ for (int i = 0; i < expected.num_rows(); ++i) {
+ for (int j = 0; j < expected.num_columns(); ++j) {
+ EXPECT_NEAR(expected_elements[i][j].real(),
+ actual_elements[i][j].real(),
+ tolerance);
+ EXPECT_NEAR(expected_elements[i][j].imag(),
+ actual_elements[i][j].imag(),
+ tolerance);
+ }
+ }
+ }
+};
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_MATRIX_TEST_HELPERS_H_
diff --git a/webrtc/modules/audio_processing/beamformer/matrix_unittest.cc b/webrtc/modules/audio_processing/beamformer/matrix_unittest.cc
new file mode 100644
index 0000000..04aed70
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/matrix_unittest.cc
@@ -0,0 +1,326 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#include <complex>
+
+#include "testing/gtest/include/gtest/gtest.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix.h"
+#include "webrtc/modules/audio_processing/beamformer/matrix_test_helpers.h"
+
+namespace webrtc {
+
+using std::complex;
+
+TEST(MatrixTest, TestMultiplySameSize) {
+ const int kNumRows = 2;
+ const int kNumCols = 2;
+ const float kValuesLeft[kNumRows][kNumCols] = {{1.1f, 2.2f}, {3.3f, 4.4f}};
+ const float kValuesRight[kNumRows][kNumCols] = {{5.4f, 127.f},
+ {4600.f, -555.f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {{10125.94f, -1081.3f},
+ {20257.82f, -2022.9f}};
+
+ Matrix<float> lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<float> rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result(kNumRows, kNumCols);
+
+ actual_result.Multiply(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEquality(expected_result, actual_result);
+
+ lh_mat.Multiply(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEquality(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestMultiplyDifferentSize) {
+ const int kNumRowsLeft = 2;
+ const int kNumColsLeft = 3;
+ const int kNumRowsRight = 3;
+ const int kNumColsRight = 2;
+ const int kValuesLeft[kNumRowsLeft][kNumColsLeft] = {{35, 466, -15},
+ {-3, 3422, 9}};
+ const int kValuesRight[kNumRowsRight][kNumColsRight] = {
+ {765, -42}, {0, 194}, {625, 66321}};
+ const int kValuesExpected[kNumRowsLeft][kNumColsRight] = {{17400, -905881},
+ {3330, 1260883}};
+
+ Matrix<int> lh_mat(*kValuesLeft, kNumRowsLeft, kNumColsLeft);
+ Matrix<int> rh_mat(*kValuesRight, kNumRowsRight, kNumColsRight);
+ Matrix<int> expected_result(*kValuesExpected, kNumRowsLeft, kNumColsRight);
+ Matrix<int> actual_result(kNumRowsLeft, kNumColsRight);
+
+ actual_result.Multiply(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEquality(expected_result, actual_result);
+
+ lh_mat.Multiply(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEquality(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestTranspose) {
+ const int kNumInitialRows = 2;
+ const int kNumInitialCols = 4;
+ const int kNumResultRows = 4;
+ const int kNumResultCols = 2;
+ const float kValuesInitial[kNumInitialRows][kNumInitialCols] = {
+ {1.1f, 2.2f, 3.3f, 4.4f}, {5.5f, 6.6f, 7.7f, 8.8f}};
+ const float kValuesExpected[kNumResultRows][kNumResultCols] = {
+ {1.1f, 5.5f}, {2.2f, 6.6f}, {3.3f, 7.7f}, {4.4f, 8.8f}};
+
+ Matrix<float> initial_mat(*kValuesInitial, kNumInitialRows, kNumInitialCols);
+ Matrix<float> expected_result(
+ *kValuesExpected, kNumResultRows, kNumResultCols);
+ Matrix<float> actual_result(kNumResultRows, kNumResultCols);
+
+ actual_result.Transpose(initial_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+ initial_mat.Transpose();
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(initial_mat, actual_result);
+}
+
+TEST(MatrixTest, TestScale) {
+ const int kNumRows = 3;
+ const int kNumCols = 3;
+ const int kScaleFactor = -9;
+ const int kValuesInitial[kNumRows][kNumCols] = {
+ {1, 20, 5000}, {-3, -29, 66}, {7654, 0, -23455}};
+ const int kValuesExpected[kNumRows][kNumCols] = {
+ {-9, -180, -45000}, {27, 261, -594}, {-68886, 0, 211095}};
+
+ Matrix<int> initial_mat(*kValuesInitial, kNumRows, kNumCols);
+ Matrix<int> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<int> actual_result;
+
+ actual_result.Scale(initial_mat, kScaleFactor);
+ MatrixTestHelpers::ValidateMatrixEquality(expected_result, actual_result);
+
+ initial_mat.Scale(kScaleFactor);
+ MatrixTestHelpers::ValidateMatrixEquality(initial_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseAdd) {
+ const int kNumRows = 2;
+ const int kNumCols = 3;
+ const float kValuesLeft[kNumRows][kNumCols] = {{1.1f, 210.45f, -549.2f},
+ {11.876f, 586.7f, -64.35f}};
+ const float kValuesRight[kNumRows][kNumCols] = {{-50.4f, 1.f, 0.5f},
+ {460.f, -554.2f, 4566.f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {
+ {-49.3f, 211.45f, -548.7f}, {471.876f, 32.5f, 4501.65f}};
+
+ Matrix<float> lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<float> rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result;
+
+ actual_result.Add(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+ lh_mat.Add(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseSubtract) {
+ const int kNumRows = 3;
+ const int kNumCols = 2;
+ const float kValuesLeft[kNumRows][kNumCols] = {
+ {1.1f, 210.45f}, {-549.2f, 11.876f}, {586.7f, -64.35f}};
+ const float kValuesRight[kNumRows][kNumCols] = {
+ {-50.4f, 1.f}, {0.5f, 460.f}, {-554.2f, 4566.f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {
+ {51.5f, 209.45f}, {-549.7f, -448.124f}, {1140.9f, -4630.35f}};
+
+ Matrix<float> lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<float> rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result;
+
+ actual_result.Subtract(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+
+ lh_mat.Subtract(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseMultiply) {
+ const int kNumRows = 1;
+ const int kNumCols = 5;
+ const float kValuesLeft[kNumRows][kNumCols] = {
+ {1.1f, 6.4f, 0.f, -1.f, -88.3f}};
+ const float kValuesRight[kNumRows][kNumCols] = {
+ {53.2f, -210.45f, -549.2f, 99.99f, -45.2f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {
+ {58.52f, -1346.88f, 0.f, -99.99f, 3991.16f}};
+
+ Matrix<float> lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<float> rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result;
+
+ actual_result.PointwiseMultiply(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+
+ lh_mat.PointwiseMultiply(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseDivide) {
+ const int kNumRows = 5;
+ const int kNumCols = 1;
+ const float kValuesLeft[kNumRows][kNumCols] = {
+ {1.1f}, {6.4f}, {0.f}, {-1.f}, {-88.3f}};
+ const float kValuesRight[kNumRows][kNumCols] = {
+ {53.2f}, {-210.45f}, {-549.2f}, {99.99f}, {-45.2f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {
+ {0.020676691f}, {-0.03041102399f}, {0.f}, {-0.010001f}, {1.9535398f}};
+
+ Matrix<float> lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<float> rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result;
+
+ actual_result.PointwiseDivide(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+
+ lh_mat.PointwiseDivide(rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(lh_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseSquareRoot) {
+ const int kNumRows = 2;
+ const int kNumCols = 2;
+ const int kValues[kNumRows][kNumCols] = {{4, 9}, {16, 0}};
+ const int kValuesExpected[kNumRows][kNumCols] = {{2, 3}, {4, 0}};
+
+ Matrix<int> operand_mat(*kValues, kNumRows, kNumCols);
+ Matrix<int> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<int> actual_result;
+
+ actual_result.PointwiseSquareRoot(operand_mat);
+ MatrixTestHelpers::ValidateMatrixEquality(expected_result, actual_result);
+
+ operand_mat.PointwiseSquareRoot();
+ MatrixTestHelpers::ValidateMatrixEquality(operand_mat, actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseSquareRootComplex) {
+ const int kNumRows = 1;
+ const int kNumCols = 3;
+ const complex<float> kValues[kNumRows][kNumCols] = {
+ {complex<float>(-4.f, 0), complex<float>(0, 9), complex<float>(3, -4)}};
+ const complex<float> kValuesExpected[kNumRows][kNumCols] = {
+ {complex<float>(0.f, 2.f), complex<float>(2.1213202f, 2.1213202f),
+ complex<float>(2.f, -1.f)}};
+
+ Matrix<complex<float> > operand_mat(*kValues, kNumRows, kNumCols);
+ Matrix<complex<float> > expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<complex<float> > actual_result;
+
+ actual_result.PointwiseSquareRoot(operand_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result,
+ actual_result);
+
+ operand_mat.PointwiseSquareRoot();
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(operand_mat,
+ actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseAbsoluteValue) {
+ const int kNumRows = 1;
+ const int kNumCols = 3;
+ const complex<float> kValues[kNumRows][kNumCols] = {
+ {complex<float>(-4.f, 0), complex<float>(0, 9), complex<float>(3, -4)}};
+ const complex<float> kValuesExpected[kNumRows][kNumCols] = {
+ {complex<float>(4.f, 0), complex<float>(9.f, 0), complex<float>(5.f, 0)}};
+
+ Matrix<complex<float> > operand_mat(*kValues, kNumRows, kNumCols);
+ Matrix<complex<float> > expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<complex<float> > actual_result;
+
+ actual_result.PointwiseAbsoluteValue(operand_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result,
+ actual_result);
+
+ operand_mat.PointwiseAbsoluteValue();
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(operand_mat,
+ actual_result);
+}
+
+TEST(MatrixTest, TestPointwiseSquare) {
+ const int kNumRows = 1;
+ const int kNumCols = 3;
+ const float kValues[kNumRows][kNumCols] = {{2.4f, -4.f, 3.3f}};
+ const float kValuesExpected[kNumRows][kNumCols] = {{5.76f, 16.f, 10.89f}};
+
+ Matrix<float> operand_mat(*kValues, kNumRows, kNumCols);
+ Matrix<float> expected_result(*kValuesExpected, kNumRows, kNumCols);
+ Matrix<float> actual_result;
+
+ actual_result.PointwiseSquare(operand_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(expected_result,
+ actual_result);
+
+ operand_mat.PointwiseSquare();
+ MatrixTestHelpers::ValidateMatrixEqualityFloat(operand_mat, actual_result);
+}
+
+TEST(MatrixTest, TestComplexOperations) {
+ const int kNumRows = 2;
+ const int kNumCols = 2;
+
+ const complex<float> kValuesLeft[kNumRows][kNumCols] = {
+ {complex<float>(1.f, 1.f), complex<float>(2.f, 2.f)},
+ {complex<float>(3.f, 3.f), complex<float>(4.f, 4.f)}};
+
+ const complex<float> kValuesRight[kNumRows][kNumCols] = {
+ {complex<float>(5.f, 5.f), complex<float>(6.f, 6.f)},
+ {complex<float>(7.f, 7.f), complex<float>(8.f, 8.f)}};
+
+ const complex<float> kValuesExpectedAdd[kNumRows][kNumCols] = {
+ {complex<float>(6.f, 6.f), complex<float>(8.f, 8.f)},
+ {complex<float>(10.f, 10.f), complex<float>(12.f, 12.f)}};
+
+ const complex<float> kValuesExpectedMultiply[kNumRows][kNumCols] = {
+ {complex<float>(0.f, 38.f), complex<float>(0.f, 44.f)},
+ {complex<float>(0.f, 86.f), complex<float>(0.f, 100.f)}};
+
+ const complex<float> kValuesExpectedPointwiseDivide[kNumRows][kNumCols] = {
+ {complex<float>(0.2f, 0.f), complex<float>(0.33333333f, 0.f)},
+ {complex<float>(0.42857143f, 0.f), complex<float>(0.5f, 0.f)}};
+
+ Matrix<complex<float> > lh_mat(*kValuesLeft, kNumRows, kNumCols);
+ Matrix<complex<float> > rh_mat(*kValuesRight, kNumRows, kNumCols);
+ Matrix<complex<float> > expected_result_add(
+ *kValuesExpectedAdd, kNumRows, kNumCols);
+ Matrix<complex<float> > expected_result_multiply(
+ *kValuesExpectedMultiply, kNumRows, kNumCols);
+ Matrix<complex<float> > expected_result_pointwise_divide(
+ *kValuesExpectedPointwiseDivide, kNumRows, kNumCols);
+ Matrix<complex<float> > actual_result_add;
+ Matrix<complex<float> > actual_result_multiply(kNumRows, kNumCols);
+ Matrix<complex<float> > actual_result_pointwise_divide;
+
+ actual_result_add.Add(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(expected_result_add,
+ actual_result_add);
+
+ actual_result_multiply.Multiply(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(
+ expected_result_multiply, actual_result_multiply);
+
+ actual_result_pointwise_divide.PointwiseDivide(lh_mat, rh_mat);
+ MatrixTestHelpers::ValidateMatrixEqualityComplexFloat(
+ expected_result_pointwise_divide, actual_result_pointwise_divide);
+}
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/pcm_utils.cc b/webrtc/modules/audio_processing/beamformer/pcm_utils.cc
new file mode 100644
index 0000000..41d729e
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/pcm_utils.cc
@@ -0,0 +1,90 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#include "webrtc/modules/audio_processing/beamformer/pcm_utils.h"
+
+#include "webrtc/base/checks.h"
+#include "webrtc/common_audio/include/audio_util.h"
+#include "webrtc/modules/audio_processing/channel_buffer.h"
+
+namespace webrtc {
+
+size_t PcmRead(FILE* file,
+ size_t length,
+ int num_channels,
+ int16_t* const* buffer) {
+ CHECK_GE(num_channels, 1);
+
+ scoped_ptr<int16_t[]> interleaved_buffer(new int16_t[length]);
+ size_t elements_read = fread(interleaved_buffer.get(), sizeof(int16_t),
+ length, file);
+ if (elements_read != length) {
+ // This is only an error if we haven't reached the end of the file.
+ CHECK_NE(0, feof(file));
+ }
+
+ Deinterleave(interleaved_buffer.get(),
+ static_cast<int>(elements_read) / num_channels,
+ num_channels,
+ buffer);
+ return elements_read;
+}
+
+size_t PcmReadToFloat(FILE* file,
+ size_t length,
+ int num_channels,
+ float* const* buffer) {
+ CHECK_GE(num_channels, 1);
+
+ int num_frames = static_cast<int>(length) / num_channels;
+ scoped_ptr<ChannelBuffer<int16_t> > deinterleaved_buffer(
+ new ChannelBuffer<int16_t>(num_frames, num_channels));
+
+ size_t elements_read =
+ PcmRead(file, length, num_channels, deinterleaved_buffer->channels());
+
+ for (int i = 0; i < num_channels; ++i) {
+ S16ToFloat(deinterleaved_buffer->channel(i), num_frames, buffer[i]);
+ }
+ return elements_read;
+}
+
+void PcmWrite(FILE* file,
+ size_t length,
+ int num_channels,
+ const int16_t* const* buffer) {
+ CHECK_GE(num_channels, 1);
+
+ scoped_ptr<int16_t[]> interleaved_buffer(new int16_t[length]);
+ Interleave(buffer,
+ static_cast<int>(length) / num_channels,
+ num_channels,
+ interleaved_buffer.get());
+ CHECK_EQ(length,
+ fwrite(interleaved_buffer.get(), sizeof(int16_t), length, file));
+}
+
+void PcmWriteFromFloat(FILE* file,
+ size_t length,
+ int num_channels,
+ const float* const* buffer) {
+ CHECK_GE(num_channels, 1);
+
+ int num_frames = static_cast<int>(length) / num_channels;
+ scoped_ptr<ChannelBuffer<int16_t> > deinterleaved_buffer(
+ new ChannelBuffer<int16_t>(num_frames, num_channels));
+
+ for (int i = 0; i < num_channels; ++i) {
+ FloatToS16(buffer[i], num_frames, deinterleaved_buffer->channel(i));
+ }
+ PcmWrite(file, length, num_channels, deinterleaved_buffer->channels());
+}
+
+} // namespace webrtc
diff --git a/webrtc/modules/audio_processing/beamformer/pcm_utils.h b/webrtc/modules/audio_processing/beamformer/pcm_utils.h
new file mode 100644
index 0000000..3a6a3b9
--- /dev/null
+++ b/webrtc/modules/audio_processing/beamformer/pcm_utils.h
@@ -0,0 +1,50 @@
+/*
+ * Copyright (c) 2014 The WebRTC project authors. All Rights Reserved.
+ *
+ * Use of this source code is governed by a BSD-style license
+ * that can be found in the LICENSE file in the root of the source
+ * tree. An additional intellectual property rights grant can be found
+ * in the file PATENTS. All contributing project authors may
+ * be found in the AUTHORS file in the root of the source tree.
+ */
+
+#ifndef WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_PCM_UTILS_H_
+#define WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_PCM_UTILS_H_
+
+#include <inttypes.h>
+#include <stdio.h>
+
+// Utilities for reading from and writing to multichannel pcm files.
+// Assumes a bit depth of 16 and little-endian. Note that in these functions,
+// length refers to the number of samples to read from/write to the file,
+// such that length / num_channels is the number of frames.
+namespace webrtc {
+
+// Reads audio from a pcm into a 2D array: buffer[channel_index][frame_index].
+// Returns the number of frames written. If this is less than |length|, it's
+// safe to assume the end-of-file was reached, as otherwise this will crash.
+// In PcmReadToFloat, the floats are within the range [-1, 1].
+size_t PcmRead(FILE* file,
+ size_t length,
+ int num_channels,
+ int16_t* const* buffer);
+size_t PcmReadToFloat(FILE* file,
+ size_t length,
+ int num_channels,
+ float* const* buffer);
+
+// Writes to a pcm file. The resulting file contains the channels interleaved.
+// Crashes if the correct number of frames aren't written to the file. For
+// PcmWriteFromFloat, floats must be within the range [-1, 1].
+void PcmWrite(FILE* file,
+ size_t length,
+ int num_channels,
+ const int16_t* const* buffer);
+void PcmWriteFromFloat(FILE* file,
+ size_t length,
+ int num_channels,
+ const float* const* buffer);
+
+} // namespace webrtc
+
+#endif // WEBRTC_MODULES_AUDIO_PROCESSING_BEAMFORMER_PCM_UTILS_H_
diff --git a/webrtc/modules/modules.gyp b/webrtc/modules/modules.gyp
index 84ec96a..3c61051 100644
--- a/webrtc/modules/modules.gyp
+++ b/webrtc/modules/modules.gyp
@@ -177,6 +177,11 @@
'audio_processing/agc/pole_zero_filter_unittest.cc',
'audio_processing/agc/standalone_vad_unittest.cc',
'audio_processing/agc/test/test_utils.cc',
+ 'audio_processing/beamformer/complex_matrix_unittest.cc',
+ 'audio_processing/beamformer/covariance_matrix_generator_unittest.cc',
+ 'audio_processing/beamformer/matrix_unittest.cc',
+ 'audio_processing/beamformer/pcm_utils.cc',
+ 'audio_processing/beamformer/pcm_utils.h',
'audio_processing/echo_cancellation_impl_unittest.cc',
'audio_processing/splitting_filter_unittest.cc',
'audio_processing/transient/dyadic_decimator_unittest.cc',