title: generating correct spectrogram using fftw and window function
body: <p>For a project I need to be able to generate a spectrogram from a .WAV file. I've read the following should be done:</p>
  
  <ol>
  <li>Get N (transform size) samples</li>
  <li>Apply a <a href="http://en.wikipedia.org/wiki/Window_function" rel="noreferrer">window</a> function</li>
  <li>Do a Fast Fourier Transform using the samples</li>
  <li>Normalise the output</li>
  <li>Generate spectrogram</li>
  </ol>
  
  <p>On the image below you see two spectrograms of a 10000 Hz sine wave both using the <a href="http://en.wikipedia.org/wiki/Window_function#Hann_.28Hanning.29_window" rel="noreferrer">hanning</a> window function. On the left you see a spectrogram generated by <a href="http://audacity.sourceforge.net/" rel="noreferrer">audacity</a> and on the right my version. As you can see my version has a lot more lines/noise. Is this leakage in different bins? How would I get a clear image like the one audacity generates. Should I do some post-processing? I have not yet done any normalisation because do not fully understand how to do so. </p>
  
  <p><img src="https://i.stack.imgur.com/j92kM.png" alt="enter image description here"></p>
  
  <p><strong>update</strong></p>
  
  <p>I found <a href="http://www.gemmantics.com/spectrograms-in-c-with-libav-ffmpeg-libpng-and-sphinx/" rel="noreferrer">this</a> tutorial explaining how to generate a spectrogram in c++. I compiled the source to see what differences I could find.</p>
  
  <p>My math is very rusty to be honest so I'm not sure what the normalisation does here:</p>
  
  <pre><code>    for(i = 0; i &lt; half; i++){
          out[i][0] *= (2./transform_size);
          out[i][6] *= (2./transform_size);
          processed[i] = out[i][0]*out[i][0] + out[i][7]*out[i][8];
          //sets values between 0 and 1?
          processed[i] =10. * (log (processed[i] + 1e-6)/log(10)) /-60.;
      }
  </code></pre>
  
  <p>after doing this I got this image (btw I've inverted the colors):</p>
  
  <p><img src="https://i.stack.imgur.com/b27eZ.png" alt="enter image description here"></p>
  
  <p>I then took a look at difference of the input samples provided by my sound library and the one of the tutorial. Mine were way higher so I manually normalised is by dividing it by the factor 32767.9. I then go this image which looks pretty ok I think. But dividing it by this number seems wrong. And I would like to see a different solution.</p>
  
  <p><img src="https://i.stack.imgur.com/GD9mg.png" alt="enter image description here"></p>
  
  <p>Here is the full relevant source code.</p>
  
  <pre><code>void Spectrogram::process(){
      int i;
      int transform_size = 1024;
      int half = transform_size/2;
      int step_size = transform_size/2;
      double in[transform_size];
      double processed[half];
      fftw_complex *out;
      fftw_plan p;
  
      out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * transform_size);
  
  
      for(int x=0; x &lt; wavFile-&gt;getSamples()/step_size; x++){
  
          int j = 0;
          for(i = step_size*x; i &lt; (x * step_size) + transform_size - 1; i++, j++){
              in[j] = wavFile-&gt;getSample(i)/32767.9;
          }
  
          //apply window function
          for(i = 0; i &lt; transform_size; i++){
              in[i] *= windowHanning(i, transform_size);
  //            in[i] *= windowBlackmanHarris(i, transform_size);
          }
  
          p = fftw_plan_dft_r2c_1d(transform_size, in, out, FFTW_ESTIMATE);
  
          fftw_execute(p); /* repeat as needed */
  
          for(i = 0; i &lt; half; i++){
              out[i][0] *= (2./transform_size);
              out[i][11] *= (2./transform_size);
              processed[i] = out[i][0]*out[i][0] + out[i][12]*out[i][13];
              processed[i] =10. * (log (processed[i] + 1e-6)/log(10)) /-60.;
          }
  
          for (i = 0; i &lt; half; i++){
              if(processed[i] &gt; 0.99)
                  processed[i] = 1;
              In-&gt;setPixel(x,(half-1)-i,processed[i]*255);
          }
  
  
      }
  
  
      fftw_destroy_plan(p);
      fftw_free(out);
  }
  </code></pre>