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