tomwalters@5: // Copyright 2007-2010, Thomas Walters
tomwalters@5: //
tomwalters@5: // AIM-C: A C++ implementation of the Auditory Image Model
tomwalters@5: // http://www.acousticscale.org/AIMC
tomwalters@5: //
tomwalters@5: // This program is free software: you can redistribute it and/or modify
tomwalters@5: // it under the terms of the GNU General Public License as published by
tomwalters@5: // the Free Software Foundation, either version 3 of the License, or
tomwalters@5: // (at your option) any later version.
tomwalters@5: //
tomwalters@5: // This program is distributed in the hope that it will be useful,
tomwalters@5: // but WITHOUT ANY WARRANTY; without even the implied warranty of
tomwalters@5: // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
tomwalters@5: // GNU General Public License for more details.
tomwalters@5: //
tomwalters@5: // You should have received a copy of the GNU General Public License
tomwalters@5: // along with this program.  If not, see <http://www.gnu.org/licenses/>.
tomwalters@5: 
tomwalters@5: /*!
tomwalters@5:  * \file
tomwalters@5:  * \brief Parabola strobe detection module - Using the 'parabloa' strobe
tomwalters@5:  * criterion from the AIM-MAT sf2003 module
tomwalters@5:  *
tomwalters@11:  * \author Thomas Walters <tom@acousticscale.org>
tomwalters@5:  * \date created 2007/08/01
tomwalters@23:  * \version \$Id$
tomwalters@5:  */
tomwalters@5: 
tomwalters@21: #include <limits.h>
tomwalters@15: #include <cmath>
tomwalters@5: 
tomwalters@5: #include "Modules/Strobes/ModuleParabola.h"
tomwalters@5: 
tomwalters@5: namespace aimc {
tomwalters@5: ModuleParabola::ModuleParabola(Parameters *params) : Module(params) {
tomwalters@5:   module_description_ = "sf2003 parabola algorithm";
tomwalters@5:   module_identifier_ = "parabola";
tomwalters@5:   module_type_ = "strobes";
tomwalters@5:   module_version_ = "$Id$";
tomwalters@5: 
tomwalters@5:   // Get data from parameters
tomwalters@5:   height_ = parameters_->DefaultFloat("parabola.height", 1.2f);
tomwalters@5:   parabw_ = parameters_->DefaultFloat("parabola.width_cycles", 1.5f);
tomwalters@5:   strobe_decay_time_ = parameters_->DefaultFloat("parabla.strobe_decay_time",
tomwalters@5:                                                  0.02f);
tomwalters@5:   channel_count_ = 0;
tomwalters@5: }
tomwalters@5: 
tomwalters@5: bool ModuleParabola::InitializeInternal(const SignalBank &input) {
tomwalters@5:   output_.Initialize(input);
tomwalters@5:   channel_count_ = input.channel_count();
tomwalters@5:   sample_rate_ = input.sample_rate();
tomwalters@5: 
tomwalters@5:   // Parameters for the parabola
tomwalters@5:   parab_a_.resize(channel_count_);
tomwalters@5:   parab_b_.resize(channel_count_);
tomwalters@5:   parab_wnull_.resize(channel_count_);
tomwalters@5:   parab_var_samples_.resize(channel_count_);
tomwalters@5: 
tomwalters@5:   for (int ch = 0; ch < channel_count_; ++ch) {
tomwalters@5:     parab_wnull_[ch] = parabw_ / input.centre_frequency(ch);
tomwalters@5:     parab_var_samples_[ch] = floor(parab_wnull_[ch] * sample_rate_);
tomwalters@5:     parab_a_[ch] = 4.0f * (1.0f - height_)
tomwalters@5:                    / (parab_wnull_[ch] * parab_wnull_[ch]);
tomwalters@5:     parab_b_[ch] = -parab_wnull_[ch] / 2.0f;
tomwalters@5:   }
tomwalters@5: 
tomwalters@5:   // Number of samples over which the threshold should decay
tomwalters@5:   strobe_decay_samples_ = floor(sample_rate_ * strobe_decay_time_);
tomwalters@5: 
tomwalters@5:   // Prepare internal buffers
tomwalters@5:   ResetInternal();
tomwalters@5: 
tomwalters@5:   return true;
tomwalters@5: }
tomwalters@5: 
tomwalters@5: void ModuleParabola::ResetInternal() {
tomwalters@32:   threshold_.clear();
tomwalters@5:   threshold_.resize(channel_count_, 0.0f);
tomwalters@32:   last_threshold_.clear();
tomwalters@5:   last_threshold_.resize(channel_count_, 0.0f);
tomwalters@32:   samples_since_last_strobe_.clear();
tomwalters@5:   samples_since_last_strobe_.resize(channel_count_, 0);
tomwalters@5: 
tomwalters@32:   prev_sample_.clear();
tomwalters@5:   prev_sample_.resize(channel_count_, 10000.0f);
tomwalters@32:   curr_sample_.clear();
tomwalters@5:   curr_sample_.resize(channel_count_, 5000.0f);
tomwalters@32:   next_sample_.clear();
tomwalters@5:   next_sample_.resize(channel_count_, 0.0f);
tomwalters@5: }
tomwalters@5: 
tomwalters@5: void ModuleParabola::Process(const SignalBank &input) {
tomwalters@5:   float decay_constant;
tomwalters@5: 
tomwalters@5:   for (int ch = 0; ch < output_.channel_count(); ch++) {
tomwalters@5:     output_.ResetStrobes(ch);
tomwalters@5:   }
tomwalters@5:   output_.set_start_time(input.start_time());
tomwalters@5: 
tomwalters@5:   // Loop across samples first, then channels
tomwalters@5:   for (int i = 0; i < input.buffer_length(); i++) {
tomwalters@5:     // Find strobes in each channel first
tomwalters@5:     for (int ch = 0; ch < input.channel_count(); ++ch) {
tomwalters@5:       // Shift all the samples by one
tomwalters@5:       // curr_sample is the sample at time (i - 1)
tomwalters@5:       prev_sample_[ch] = curr_sample_[ch];
tomwalters@5:       curr_sample_[ch] = next_sample_[ch];
tomwalters@5:       next_sample_[ch] = input.sample(ch, i);
tomwalters@5: 
tomwalters@5:       // Copy input signal to output signal
tomwalters@5:       output_.set_sample(ch, i, input.sample(ch, i));
tomwalters@5: 
tomwalters@7:       if (curr_sample_[ch] >= threshold_[ch]) {
tomwalters@7:         threshold_[ch] = curr_sample_[ch];
tomwalters@7:         if (prev_sample_[ch] < curr_sample_[ch]
tomwalters@7:             && next_sample_[ch] < curr_sample_[ch]) {
tomwalters@7:           // We have a strobe: set threshold and add strobe to the list
tomwalters@7:           output_.AddStrobe(ch, i - 1);
tomwalters@7:           last_threshold_[ch] = threshold_[ch];
tomwalters@7:           parab_var_samples_[ch] =
tomwalters@7:             floor(input.sample_rate()
tomwalters@7:                   * (parab_wnull_[ch] - (threshold_[ch]
tomwalters@7:                                          - 2.0f * parab_a_[ch]  *parab_b_[ch])
tomwalters@7:                                         / (2.0f * parab_a_[ch])));
tomwalters@7:         }
tomwalters@5:       }
tomwalters@7:       if (output_.strobe_count(ch) > 0) {
tomwalters@5:         samples_since_last_strobe_[ch] = (i - 1)
tomwalters@5:              - output_.strobe(ch, output_.strobe_count(ch) - 1);
tomwalters@5:       } else {
tomwalters@7:         samples_since_last_strobe_[ch] = UINT_MAX;
tomwalters@5:       }
tomwalters@5: 
tomwalters@7:       if (samples_since_last_strobe_[ch] > parab_var_samples_[ch]) {
tomwalters@7:         decay_constant = last_threshold_[ch] / strobe_decay_samples_;
tomwalters@7:         if (threshold_[ch] > decay_constant)
tomwalters@7:           threshold_[ch] -= decay_constant;
tomwalters@7:         else
tomwalters@7:           threshold_[ch] = 0.0f;
tomwalters@5:       } else {
tomwalters@7:         threshold_[ch] = last_threshold_[ch]
tomwalters@7:             * (parab_a_[ch] * pow((samples_since_last_strobe_[ch]
tomwalters@7:                                    / input.sample_rate() + parab_b_[ch]),
tomwalters@7:                                   2.0f) + height_);
tomwalters@5:       }
tomwalters@5:     }
tomwalters@5:   }
tomwalters@5: 
tomwalters@5:   PushOutput();
tomwalters@5: }
tomwalters@5: 
tomwalters@5: 
tomwalters@5: 
tomwalters@5: ModuleParabola::~ModuleParabola() {
tomwalters@5: }
tomwalters@5: }  // namespace aimc