/* Program untuk membandingkan besar nois antara downsamplig dengan 
   antialiasing filter dengan dan tanpa antialiasing filter */
/* hasil dari program ini adalah aliasing.raw dan antialiasing.raw */
/* silahkan bandingkan dengan cara play kedua file tersebut */

#include <stdio.h> 
#include <stdlib.h>
#include <string.h> 
#include <alsa/asoundlib.h> 
#include <unistd.h> 
#include <math.h> 
#include <time.h> 


int main(int argc, char *argv[])
{
  int i;
  char *file_name;
  file_name = malloc(sizeof(char) * 50); 
  int rtime =114 ; // 114 detik estimasi jitter transmisi WSPR
  int buffer_frames_48 = rtime*48000; // alokasi buffer sampling 48KHz sesuai jumlah sampel
  int buffer_frames_12 = rtime*12000; // alokasi buffer sampleing 12 KHz
  //  unsigned int rate = 48000;// sampling rate
  int16_t *buf2_int16, *buf2;   // buffer integer 16bits
  float *buf2_float;


  if(argc == 1){
    printf("\nkurang namafile *.raw...\n");
    exit(-1);
  }
    if (argc > 1)
    file_name = argv[1];
    printf("\n%s\n",argv[1]);

    buf2 = malloc(buffer_frames_48 * 4);
    buf2_int16 = malloc(buffer_frames_48 * 2);
    buf2_float = malloc(buffer_frames_48*sizeof(float));


    FILE *fin = fopen(file_name, "rb"); if (fin == NULL) {
      printf("Kesalahan membuka file!\n");
      exit(-1);
    }

    if((fread(buf2, 4, buffer_frames_48, fin)) == 0)exit(-1);       /* jika stereo */
    fclose(fin); 

    /* ambil 1 sampel mono  */

    for (i = 0; i < buffer_frames_48; i++){
      buf2[i]=buf2[2*i];
    }

    /**** decimasi menjadi 12 KHz sampling ***/

    for (i = 0; i < buffer_frames_12; i++){
      buf2_int16[i] = buf2[i*4];
    }

    double pwr = 0;
    for (i = 0; i < buffer_frames_12; i++){
      pwr=pwr+ (double)(buf2_int16[i]*buf2_int16[i]);
    }
    pwr=pwr/(double)buffer_frames_12;
    pwr=10*log10(pwr);
    printf("Level rerata power aliasing = %.1lf dB\n",pwr);
    
    FILE *fal = fopen("./alias.WAV", "wb"); if (fal == NULL) {
      printf("Kesalahan membuka file!\n");
      exit(-11);
    }
    

    short header[22];
    for(i=0;i<22;i++)header[i]=0;
    fwrite(header, 2, 22, fal);
    fwrite(buf2_int16, 2, buffer_frames_12, fal);
    fclose(fal); 

    /* anti aliasing filter  */
    
    float pi = 3.141592654; 
    float fc = 2500;
    float frq=48000;
    float a1,b1,c1,tau1;

  tau1=1/(2*pi*fc);
  a1=(1.0-1.0/(2.0*tau1*frq))/(1.0+1.0/(2.0*tau1*frq));
  b1=c1=(1.0/(2.0*tau1*frq))/(1.0+1.0/(2.0*tau1*frq));

  buf2_int16[0]=buf2[0];
  buf2_float[0]=(float)buf2[0];
  for (i = 1; i < buffer_frames_48; i++){
    buf2_float[i]=a1*buf2_float[i-1]+b1*(float)buf2[i]+c1*(float)buf2[i-1];
    buf2_int16[i]=(int16_t)buf2_float[i];
  }
  /*************************/

  /**** decimasi menjadi 12 KHz sampling ***/

  for (i = 0; i < buffer_frames_12; i++){
    buf2_int16[i] = buf2_int16[i*4];
  }

    pwr = 0;
    for (i = 0; i < buffer_frames_12; i++){
      pwr=pwr+ (double)(buf2_int16[i]*buf2_int16[i]);
    }
    pwr=pwr/(double)buffer_frames_12;
    pwr=10*log10(pwr);
    printf("Level rerata power antialiasing = %.1lf dB\n",pwr);

    FILE *fan = fopen("./antialias.WAV", "wb"); if (fan == NULL) {
    printf("Kesalahan membuka file!\n");
    exit(1);
  }

  fwrite(header, 2, 22, fan);
  fwrite(buf2_int16, 2, buffer_frames_12, fan);
  fclose(fan); 


  free(buf2);
  free(buf2_int16);
  free(buf2_float);
}