付録D 相関関数を求めるプログラム

/* 周波数に変換したもので相関関数を求める */

#include        <stdio.h>
#include        <stdlib.h>
#include        <string.h>
#include        <math.h>

#define   MAX   1000
#define   filenumber    1

int  tsmax1,tsmax2;
double  base[MAX],arange[MAX];
FILE  *fpin1,*fpin2,*fpout;

void sougosoukan();

void main( int argc , char *argv[] ){

  int  i,j;
  char  infile1[MAX],infile2[MAX];

  if( argc != 3 ){
    fprintf( stderr,"usage: soukan infile infile2\n" );
    exit(1);
  }

  strcpy( infile1 , argv[1] );
  strcpy( infile2 , argv[2] );

  if(( fpin1 = fopen( infile1,"r" )) == NULL ){
    fprintf( stderr,"Can't Open infile1 \n" );
    exit(1);
  }
  if(( fpin2 = fopen( infile2,"r" )) == NULL ){
    fprintf( stderr,"Can't Open infile2 \n" );
    exit(1);
  }

  /* 最大時間数を求めておく */
  i = 0;
  while( fscanf( fpin1 , "%lf",&base[i] ) == 1 ){
    i++;
  }
  tsmax1 = i;
  j = 0;
  while( fscanf( fpin2 , "%lf",&arange[j] ) == 1 ){
    j++;
  }
  tsmax2 = j;

  sougosoukan();

  printf( "tsmax1 =%d\n",tsmax1 );
  printf( "tsmax2 =%d\n",tsmax2 );

  fclose( fpin1 );
  fclose( fpin2 );
  
}

void   sougosoukan(){

  int        i,t,l;
  int        length;
  double     sumbase,sumarange,sumsougo;
  double     cbase[MAX],carange[MAX],csougo[MAX];
  char  outfile[MAX];

  sprintf( outfile ,"syusoukan_ver%d.dat" ,filenumber );
  fpout = fopen( outfile ,"w" );
  if( fpout == NULL ){
    fprintf( stderr, "File(output) not open!\n" );
    exit(1);
  }

  for( t = 0 ; t < tsmax1 ; t++ ){
    base[t + tsmax1] = base[t];
  }
  for( t = 0 ; t < tsmax2 ; t++ ){
    arange[t + tsmax2] = arange[t];
  }

  /* 原曲の自己相関 */
  for( l = 0 ; l < tsmax1 ; l++ ){
    sumbase = 0;
    for( i = 0 ; i < tsmax1 ; i++ ){
      sumbase += base[i] * base[i + l];
    }
    cbase[l] = sumbase / tsmax1;
  }

  /* 編曲された曲の自己相関 */
  for( l = 0 ; l < tsmax2 ; l++ ){
    sumarange = 0;
    for( i = 0 ; i < tsmax2 ; i++ ){
      sumarange += arange[i] * arange[i + l];
    }
    carange[l] = sumarange / tsmax2;
  }

  /* ２つの曲の相互相関 */
  if( tsmax1 < tsmax2 ){
    length = tsmax2;
  }else{
    length = tsmax1;
  }
  for( l = 0 ; l < length ; l++ ){
    sumsougo = 0;
    for( i = 0 ; i < length ; i++ ){
      sumsougo += base[i] * arange[i + l];
    }
    csougo[l] = sumsougo / length;
  }

  /* ファイルに出力する */
  for( t = 0 ; t < length ; t++ ){
    fprintf( fpout , "base[%d] = %f",t,cbase[t] );
    fprintf( fpout , " arange[%d] = %f",t,carange[t] );
    fprintf( fpout , " sougo[%d] = %f\n",t,csougo[t] );
  }

  fclose( fpout );
}