shun_iwasawa a35b8f
/*
shun_iwasawa a35b8f
Copyright (c) 2003-2004, Mark Borgerding
shun_iwasawa a35b8f
shun_iwasawa a35b8f
All rights reserved.
shun_iwasawa a35b8f
shun_iwasawa a35b8f
Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    * Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
shun_iwasawa a35b8f
    * Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
shun_iwasawa a35b8f
    * Neither the author nor the names of any contributors may be used to endorse or promote products derived from this software without specific prior written permission.
shun_iwasawa a35b8f
shun_iwasawa a35b8f
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
shun_iwasawa a35b8f
*/
shun_iwasawa a35b8f
shun_iwasawa a35b8f
#include <stdlib.h></stdlib.h>
shun_iwasawa a35b8f
#include <math.h></math.h>
shun_iwasawa a35b8f
#include <stdio.h></stdio.h>
shun_iwasawa a35b8f
#include <string.h></string.h>
shun_iwasawa a35b8f
#include <unistd.h></unistd.h>
shun_iwasawa a35b8f
#include <png.h></png.h>
shun_iwasawa a35b8f
shun_iwasawa a35b8f
#include "kiss_fft.h"
shun_iwasawa a35b8f
#include "kiss_fftr.h"
shun_iwasawa a35b8f
shun_iwasawa a35b8f
int nfft=1024;
shun_iwasawa a35b8f
FILE * fin=NULL;
shun_iwasawa a35b8f
FILE * fout=NULL;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
int navg=20;
shun_iwasawa a35b8f
int remove_dc=0;
shun_iwasawa a35b8f
int nrows=0;
shun_iwasawa a35b8f
float * vals=NULL;
shun_iwasawa a35b8f
int stereo=0;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
static
shun_iwasawa a35b8f
void config(int argc,char** argv)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    while (1) {
shun_iwasawa a35b8f
        int c = getopt (argc, argv, "n:r:as");
shun_iwasawa a35b8f
        if (c == -1)
shun_iwasawa a35b8f
            break;
shun_iwasawa a35b8f
        switch (c) {
shun_iwasawa a35b8f
        case 'n': nfft=(int)atoi(optarg);break;
shun_iwasawa a35b8f
        case 'r': navg=(int)atoi(optarg);break;
shun_iwasawa a35b8f
        case 'a': remove_dc=1;break;
shun_iwasawa a35b8f
        case 's': stereo=1;break;
shun_iwasawa a35b8f
        case '?':
shun_iwasawa a35b8f
            fprintf (stderr, "usage options:\n"
shun_iwasawa a35b8f
                     "\t-n d: fft dimension(s) [1024]\n"
shun_iwasawa a35b8f
                     "\t-r d: number of rows to average [20]\n"
shun_iwasawa a35b8f
                     "\t-a : remove average from each fft buffer\n"
shun_iwasawa a35b8f
                     "\t-s : input is stereo, channels will be combined before fft\n"
shun_iwasawa a35b8f
                     "16 bit machine format real input is assumed\n"
shun_iwasawa a35b8f
                     );
shun_iwasawa a35b8f
        default:
shun_iwasawa a35b8f
            fprintf (stderr, "bad %c\n", c);
shun_iwasawa a35b8f
            exit (1);
shun_iwasawa a35b8f
            break;
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
    if ( optind < argc ) {
shun_iwasawa a35b8f
        if (strcmp("-",argv[optind]) !=0)
shun_iwasawa a35b8f
            fin = fopen(argv[optind],"rb");
shun_iwasawa a35b8f
        ++optind;
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    if ( optind < argc ) {
shun_iwasawa a35b8f
        if ( strcmp("-",argv[optind]) !=0 ) 
shun_iwasawa a35b8f
            fout = fopen(argv[optind],"wb");
shun_iwasawa a35b8f
        ++optind;
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
    if (fin==NULL)
shun_iwasawa a35b8f
        fin=stdin;
shun_iwasawa a35b8f
    if (fout==NULL)
shun_iwasawa a35b8f
        fout=stdout;
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
#define CHECKNULL(p) if ( (p)==NULL ) do { fprintf(stderr,"CHECKNULL failed @ %s(%d): %s\n",__FILE__,__LINE__,#p );exit(1);} while(0)
shun_iwasawa a35b8f
shun_iwasawa a35b8f
typedef struct
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    png_byte r;
shun_iwasawa a35b8f
    png_byte g;
shun_iwasawa a35b8f
    png_byte b;
shun_iwasawa a35b8f
} rgb_t;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
static 
shun_iwasawa a35b8f
void val2rgb(float x,rgb_t *p)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    const double pi = 3.14159265358979;
shun_iwasawa a35b8f
    p->g = (int)(255*sin(x*pi));
shun_iwasawa a35b8f
    p->r = (int)(255*abs(sin(x*pi*3/2)));
shun_iwasawa a35b8f
    p->b = (int)(255*abs(sin(x*pi*5/2)));
shun_iwasawa a35b8f
    //fprintf(stderr,"%.2f : %d,%d,%d\n",x,(int)p->r,(int)p->g,(int)p->b);
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
static
shun_iwasawa a35b8f
void cpx2pixels(rgb_t * res,const float * fbuf,size_t n)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    size_t i;
shun_iwasawa a35b8f
    float minval,maxval,valrange;
shun_iwasawa a35b8f
    minval=maxval=fbuf[0];
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    for (i = 0; i < n; ++i) {
shun_iwasawa a35b8f
        if (fbuf[i] > maxval) maxval = fbuf[i];
shun_iwasawa a35b8f
        if (fbuf[i] < minval) minval = fbuf[i];
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    fprintf(stderr,"min ==%f,max=%f\n",minval,maxval);
shun_iwasawa a35b8f
    valrange = maxval-minval;
shun_iwasawa a35b8f
    if (valrange == 0) {
shun_iwasawa a35b8f
        fprintf(stderr,"min == max == %f\n",minval);
shun_iwasawa a35b8f
        exit (1);
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    for (i = 0; i < n; ++i)
shun_iwasawa a35b8f
        val2rgb( (fbuf[i] - minval)/valrange , res+i );
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
static
shun_iwasawa a35b8f
void transform_signal(void)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    short *inbuf;
shun_iwasawa a35b8f
    kiss_fftr_cfg cfg=NULL;
shun_iwasawa a35b8f
    kiss_fft_scalar *tbuf;
shun_iwasawa a35b8f
    kiss_fft_cpx *fbuf;
shun_iwasawa a35b8f
    float *mag2buf;
shun_iwasawa a35b8f
    int i;
shun_iwasawa a35b8f
    int n;
shun_iwasawa a35b8f
    int avgctr=0;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    int nfreqs=nfft/2+1;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    CHECKNULL( cfg=kiss_fftr_alloc(nfft,0,0,0) );
shun_iwasawa a35b8f
    CHECKNULL( inbuf=(short*)malloc(sizeof(short)*2*nfft ) );
shun_iwasawa a35b8f
    CHECKNULL( tbuf=(kiss_fft_scalar*)malloc(sizeof(kiss_fft_scalar)*nfft ) );
shun_iwasawa a35b8f
    CHECKNULL( fbuf=(kiss_fft_cpx*)malloc(sizeof(kiss_fft_cpx)*nfreqs ) );
shun_iwasawa a35b8f
    CHECKNULL( mag2buf=(float*)malloc(sizeof(float)*nfreqs ) );
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    memset(mag2buf,0,sizeof(mag2buf)*nfreqs);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    while (1) {
shun_iwasawa a35b8f
        if (stereo) {
shun_iwasawa a35b8f
            n = fread(inbuf,sizeof(short)*2,nfft,fin);
shun_iwasawa a35b8f
            if (n != nfft ) 
shun_iwasawa a35b8f
                break;
shun_iwasawa a35b8f
            for (i=0;i
shun_iwasawa a35b8f
                tbuf[i] = inbuf[2*i] + inbuf[2*i+1];
shun_iwasawa a35b8f
        }else{
shun_iwasawa a35b8f
            n = fread(inbuf,sizeof(short),nfft,fin);
shun_iwasawa a35b8f
            if (n != nfft ) 
shun_iwasawa a35b8f
                break;
shun_iwasawa a35b8f
            for (i=0;i
shun_iwasawa a35b8f
                tbuf[i] = inbuf[i];
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
        if (remove_dc) {
shun_iwasawa a35b8f
            float avg = 0;
shun_iwasawa a35b8f
            for (i=0;i
shun_iwasawa a35b8f
            avg /= nfft;
shun_iwasawa a35b8f
            for (i=0;i
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
        /* do FFT */
shun_iwasawa a35b8f
        kiss_fftr(cfg,tbuf,fbuf);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
        for (i=0;i
shun_iwasawa a35b8f
            mag2buf[i] += fbuf[i].r * fbuf[i].r + fbuf[i].i * fbuf[i].i;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
        if (++avgctr == navg) {
shun_iwasawa a35b8f
            avgctr=0;
shun_iwasawa a35b8f
            ++nrows;
shun_iwasawa a35b8f
            vals = (float*)realloc(vals,sizeof(float)*nrows*nfreqs);
shun_iwasawa a35b8f
            float eps = 1;
shun_iwasawa a35b8f
            for (i=0;i
shun_iwasawa a35b8f
                vals[(nrows - 1) * nfreqs + i] = 10 * log10 ( mag2buf[i] / navg + eps );
shun_iwasawa a35b8f
            memset(mag2buf,0,sizeof(mag2buf[0])*nfreqs);
shun_iwasawa a35b8f
        }
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    free(cfg);
shun_iwasawa a35b8f
    free(inbuf);
shun_iwasawa a35b8f
    free(tbuf);
shun_iwasawa a35b8f
    free(fbuf);
shun_iwasawa a35b8f
    free(mag2buf);
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
static
shun_iwasawa a35b8f
void make_png(void)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    png_bytepp row_pointers=NULL;
shun_iwasawa a35b8f
    rgb_t * row_data=NULL;
shun_iwasawa a35b8f
    int i;
shun_iwasawa a35b8f
    int nfreqs = nfft/2+1;
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    png_structp png_ptr=NULL;
shun_iwasawa a35b8f
    png_infop info_ptr=NULL;
shun_iwasawa a35b8f
    
shun_iwasawa a35b8f
    CHECKNULL( png_ptr = png_create_write_struct (PNG_LIBPNG_VER_STRING,0,0,0) );
shun_iwasawa a35b8f
    CHECKNULL( info_ptr = png_create_info_struct(png_ptr) );
shun_iwasawa a35b8f
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    png_init_io(png_ptr, fout );
shun_iwasawa a35b8f
    png_set_IHDR(png_ptr, info_ptr ,nfreqs,nrows,8,PNG_COLOR_TYPE_RGB,PNG_INTERLACE_NONE,PNG_COMPRESSION_TYPE_DEFAULT,PNG_FILTER_TYPE_DEFAULT );
shun_iwasawa a35b8f
    
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    row_data = (rgb_t*)malloc(sizeof(rgb_t) * nrows * nfreqs) ;
shun_iwasawa a35b8f
    cpx2pixels(row_data, vals, nfreqs*nrows );
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    row_pointers = realloc(row_pointers, nrows*sizeof(png_bytep));
shun_iwasawa a35b8f
    for (i=0;i
shun_iwasawa a35b8f
        row_pointers[i] = (png_bytep)(row_data + i*nfreqs);
shun_iwasawa a35b8f
    }
shun_iwasawa a35b8f
    png_set_rows(png_ptr, info_ptr, row_pointers);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    fprintf(stderr,"creating %dx%d png\n",nfreqs,nrows);
shun_iwasawa a35b8f
    fprintf(stderr,"bitdepth %d \n",png_get_bit_depth(png_ptr,info_ptr ) );
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    png_write_png(png_ptr, info_ptr, PNG_TRANSFORM_IDENTITY , NULL);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
}
shun_iwasawa a35b8f
shun_iwasawa a35b8f
int main(int argc,char ** argv)
shun_iwasawa a35b8f
{
shun_iwasawa a35b8f
    config(argc,argv);
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    transform_signal();
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    make_png();
shun_iwasawa a35b8f
shun_iwasawa a35b8f
    if (fout!=stdout) fclose(fout);
shun_iwasawa a35b8f
    if (fin!=stdin) fclose(fin);
shun_iwasawa a35b8f
    return 0;
shun_iwasawa a35b8f
}