Bandlimited square-wave synthesis using additive techniques.

bandlimitedSquareWaveGeneration.cxx

/*
 * Author: Harry van Haaren 2013
 *         harryhaaren@gmail.com
 * 
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
 * MA 02110-1301, USA.
 */

#include <math.h>
#include <sstream>
#include <fstream>
#include <iostream>

#include <fstream>
#include <string>
#include <iostream>

#define SAMPLE_RATE 44100

void square(int count, float* array, int freq, int harmonics)
{
  for(float i = 0.f; i < count; i++ )
  {
    // calculate wavelength (in samples) of frequency
    float wavelength = SAMPLE_RATE / freq;
    
    // temporary value, for adding harmonics before writing to array
    float out = 0.f;
    
    // loop over each harmonic
    for(int h = 1; h < harmonics + 1; h++ )
    {
      float harmonicFreq = (2*h-1) * freq;
      
      // check harmonic for aliasing ( > nyquist)
      if ( harmonicFreq < SAMPLE_RATE / 2.f )
      {
        // formula for square wave:
        //    odd harmonics only
        //    harmonic amplitude  = 1 / <harmonic number>
        out += sin( (2*h-1) * (i/wavelength) ) / (2*h-1);
      }
    }
    
    // write the output value of the added harmonics to the array
    array[int(i)] = out / 2.f;
  }
}

int main()
{
  
  int freqs[] = {110, 220, 440, 880, 1760, 3520, 7040, 14080};
  int harms[8];
  
  std::stringstream output;
  output << "namespace square {";
  
  for( int i = 0; i < 8; i++)
  {
    harms[i] = SAMPLE_RATE / freqs[i];
    printf("f %i = h %i\n", freqs[i], harms[i] );
    
    int wavetableSize = SAMPLE_RATE / freqs[i] + 1;
    
    // setup array
    float array[wavetableSize];
    
    // do the calculation
    square( wavetableSize, &array[0], freqs[i], harms[i] );
    
    // write the array to file
    std::stringstream s;
    s << "oct" << i;
    
    output << "float " << s.str() << "[] = {";
    
    for(int i = 0; i < wavetableSize; i++)
    {
      if( i % 10 == 0 )
      {
        output << "\n";
      }
      
      output << array[i] << ", ";
      
    }
    
    output << "};" ;
  }
  
  output << "\n}\n\n";
  
  std::ofstream write;
  write.open( "square.hxx" );
  write << output.str();
  write.close();
}