sofia@0: /*
sofia@0:   GPU multi-rate FIR filter bank example software
sofia@0: 
sofia@0:   Oxford e-Research Centre, Oxford University
sofia@0: 
sofia@0:   Centre for Digital Music, Queen Mary, University of London.
sofia@0: 
sofia@0:   This program is free software: you can redistribute it and/or modify
sofia@0:   it under the terms of the GNU General Public License as published by
sofia@0:   the Free Software Foundation, either version 3 of the License, 
sofia@0:   or (at your option) any later version.
sofia@0:   See the file COPYING included with this distribution for more information.
sofia@0: */
sofia@0: 
sofia@0: #include <cstdio>
sofia@0: #include <cstdlib>
sofia@0: #include <cstring>
sofia@0: 
sofia@0: #include <sys/time.h>
sofia@0: #include "filters.h"
sofia@0: 
sofia@0: ////////////////////////////////////////////////////////////////////////////////
sofia@0: // Program main
sofia@0: ////////////////////////////////////////////////////////////////////////////////
sofia@0: 
sofia@0: using namespace std;
sofia@0: 
sofia@0: 
sofia@0: int main( int argc, char** argv) 
sofia@0: {
sofia@0:   unsigned int N = 100; // total number of input blocks before resampling
sofia@0: 
sofia@0:   int bufb_size = B_SIZE + OFFSET;
sofia@0: 
sofia@0: // the total input signal to loop through is
sofia@0: // essentially used to check results and provide multiple
sofia@0: // block iterations.
sofia@0: 
sofia@0:   float *h_in[MAX_RATES], *h_out[MAX_RATES], *h_reference[MAX_RATES]; // input and output arrays
sofia@0:   int out_sz[MAX_RATES];  //total output size for all filters
sofia@0: 
sofia@0:   //timing
sofia@0:   struct timeval t_start, t_end;
sofia@0:   double t_gpu = 0.0, t_mcpu = 0.0;
sofia@0: 
sofia@0:   // parameters, filter and GPU arrays structures
sofia@0:   // object pointers can be passed to these for C++ code
sofia@0:   filter_arrays farr;
sofia@0:   params gparams;
sofia@0:   gpu_arrays gpuarrays;
sofia@0:   int ratenumf[MAX_RATES];
sofia@0:  
sofia@0:   // command line arguments
sofia@0:   cmd_args args;
sofia@0: 
sofia@0:   //deaults
sofia@0:   args.nrates = 3;
sofia@0:   args.nf = 60;
sofia@0:   args.insize = 1024;
sofia@0:   args.rconst = 0;
sofia@0:   args.tim=1;
sofia@0: 
sofia@0:   read_command_line(argc, argv, &args);
sofia@0: 
sofia@0:   //initialise parameters
sofia@0:   int nf = args.nf;
sofia@0:   int numrates = args.nrates;
sofia@0:   int rem= nf%numrates;
sofia@0: 
sofia@0:   printf("\nGPU FIR filter parameters\n------\n");
sofia@0:   printf("\nTotal number of input blocks = %d\n", N);
sofia@0:   gparams.nfilters = nf; // not sure if this is needed
sofia@0:   printf("\nTotal number of filters = %d\n", gparams.nfilters);
sofia@0:   gparams.fsize = B_SIZE;
sofia@0:   printf("\nFilter size = %d\n", gparams.fsize);
sofia@0:   gparams.nrates = numrates;
sofia@0:   printf("\nNumber of sampling rates = %d\n", gparams.nrates);
sofia@0:   gparams.streams = 1;
sofia@0:   printf("\nCUDA streams flag = %d\n", gparams.streams);
sofia@0: 
sofia@0:   //dividing sampling rates equally...
sofia@0:   for (int i=0; i< numrates; ++i){
sofia@0:    ratenumf[i]  = nf/numrates;
sofia@0:   }
sofia@0:   if (rem > 0){
sofia@0:     for (int i=0; i<rem; ++i)
sofia@0:       ratenumf[i]++;
sofia@0:   }
sofia@0: 
sofia@0:   int in_sz = args.insize; 
sofia@0:   if (numrates == 1 && args.rconst==1){
sofia@0:       gpuarrays.osize[0] = in_sz;
sofia@0:   }
sofia@0:   else{
sofia@0:     for (int i=0; i< numrates; ++i){
sofia@0:       gpuarrays.osize[i] = in_sz/rdiv[i];
sofia@0:       if (gpuarrays.osize[i] < B_SIZE ){
sofia@0: 	printf("\n\nInput size for rate %d is  shorter than filter size.\nChose a longer input block or a shorter filter.\n", i);
sofia@0: 	printf("\nFilter size = %d\n", B_SIZE);
sofia@0: 	printf("\nDecimated  size = %d\n", gpuarrays.osize[i]);
sofia@0: 	printf("\nInitial input block size = %d\n", in_sz);
sofia@0: 	printf("\nDecimation factor = %d\n\n", rdiv[i]);
sofia@0: 	exit(EXIT_FAILURE);
sofia@0:       }
sofia@0:     }
sofia@0:   }
sofia@0: 
sofia@0:   int out_blk_sz = 0;
sofia@0: 
sofia@0:   for (int i=0; i<numrates; i++){
sofia@0:     gparams.rnumf[i] = ratenumf[i];
sofia@0:     printf("\nNumber of filters for rate %d = %d\n", i, gparams.rnumf[i]);
sofia@0:     gpuarrays.isize[i] = gpuarrays.osize[i] + gparams.fsize;
sofia@0:     printf("\nGPU Input size for rate %d = %d\n", i, gpuarrays.isize[i]);
sofia@0:     out_sz[i] = N * gpuarrays.osize[i];
sofia@0:     out_blk_sz +=ratenumf[i] *gpuarrays.osize[i];
sofia@0:     printf("\nGPU output size for 1 block for this rate: %d\n", ratenumf[i] *gpuarrays.osize[i]);
sofia@0:   }
sofia@0:     printf("\ntotal output size for 1 block for all rates: %d\n", out_blk_sz);
sofia@0:     printf("----------------\n");
sofia@0: 
sofia@0:   int oindex = 0;
sofia@0:   int pos[MAX_RATES];
sofia@0:   int sum = 0;
sofia@0:   pos[0] = 0;
sofia@0:   for (int r=1; r<numrates; r++){
sofia@0:     sum+=ratenumf[r-1];
sofia@0:     pos[r] = sum;
sofia@0:   }
sofia@0: 
sofia@0:   // Initialize arrays
sofia@0:   for (int i=0; i<numrates; i++){
sofia@0:     h_in[i] = (float*) malloc((gpuarrays.osize[i]*N + gparams.fsize)*sizeof(float));
sofia@0: 
sofia@0:     h_out[i] = (float*) malloc(gpuarrays.osize[i]*ratenumf[i]*N*sizeof(float));
sofia@0:     
sofia@0:     h_reference[i] = (float*) malloc(gpuarrays.osize[i]*ratenumf[i]*N*sizeof(float));
sofia@0: 
sofia@0:     for (int n=0; n<gparams.fsize; ++n)
sofia@0:       h_in[i][n] = 0.0f;
sofia@0: 
sofia@0:     for (int n=0; n<gpuarrays.osize[i]*N; ++n)
sofia@0:       h_in[i][n + gparams.fsize] = rand() / (float)RAND_MAX; 
sofia@0: 
sofia@0:     for (int n=0; n<gpuarrays.osize[i]*ratenumf[i]*N; ++n){
sofia@0:       h_out[i][n] = 0.0f;
sofia@0:       h_reference[i][n] = 0.0f;
sofia@0:     }
sofia@0:   }
sofia@0: 
sofia@0:   // initialize filters with random numbers
sofia@0:   farr.bk = (float*) malloc(nf*gparams.fsize*sizeof(float)); 
sofia@0:   for (int i=0; i < nf; ++i){
sofia@0:     farr.m_offb[i] = OFFSET;
sofia@0:     for (int b=0; b < gparams.fsize;  ++b)
sofia@0:       farr.bk[i*gparams.fsize+b] = rand() / (float)RAND_MAX;
sofia@0:   }
sofia@0:   
sofia@0:   //initialize cpu buffers for each filter 
sofia@0:   for (int f=0; f < nf; ++f){
sofia@0:     for (int i=0; i < bufb_size ; i++){
sofia@0:       farr.buf_in[f][i] = 0.0;
sofia@0:     }
sofia@0:   }	  
sofia@0: 
sofia@0:   printf("\nRunning multirate filter bank test on GPU, CPU, and CPU OpenMP\n");
sofia@0:   //  printf("----------------\n");
sofia@0:   // compute reference solution  
sofia@0:   printf("\ncompute reference solution (1 CPU)\n");
sofia@0: 
sofia@0:   compute_ref( h_in, h_reference, &gpuarrays, &gparams, &args, &farr, N);
sofia@0: 
sofia@0:   //reset CPU buffers
sofia@0:   for (int f=0; f < nf; ++f){
sofia@0:     farr.m_offb[f] = OFFSET;
sofia@0:     for (int i=0; i < bufb_size ; ++i)
sofia@0:       farr.buf_in[f][i] = 0.0;
sofia@0:   }
sofia@0: 
sofia@0: // compute OPENMP solution    
sofia@0:   printf("\nCompute OPENMP solution\n");
sofia@0: 
sofia@0:   compute_omp( h_in, h_reference, &gpuarrays, &gparams, &args, &farr, N);
sofia@0: 
sofia@0:  // compute CUDA solution  
sofia@0:  printf("\nCompute cuda solution\n");
sofia@0:  
sofia@0: // init
sofia@0:  cudaMultiFilterFirInit(h_in, h_out, farr.bk, &gparams, &gpuarrays);
sofia@0:  for (int i=0; i<numrates; i++)
sofia@0:    memcpy(gpuarrays.h_in[i], &h_in[i][0], gpuarrays.isize[i]*sizeof(float) );
sofia@0:  //time execution
sofia@0:  oindex = 0;
sofia@0:  if (args.tim)
sofia@0:    gettimeofday(&t_start, NULL);
sofia@0:  
sofia@0:  for (int ii=0; ii<N; ++ii){ // loop through input blocks
sofia@0:    
sofia@0:    //need to copy data to pinned memory for streams...
sofia@0:    for (int i=0; i<numrates; i++){
sofia@0:      if (ii==0)
sofia@0:        memcpy(gpuarrays.h_in[i], &h_in[i][0], gpuarrays.isize[i]*sizeof(float) );
sofia@0:      else
sofia@0:        memcpy(gpuarrays.h_in[i], &h_in[i][ii*gpuarrays.osize[i]], gpuarrays.isize[i]*sizeof(float) );
sofia@0: }
sofia@0: 
sofia@0:    //call GPU function   
sofia@0:    cudaMultiFilterFirStreams(&gpuarrays, &gparams);
sofia@0:    
sofia@0:    // ... and copy data back from pinned memory...
sofia@0:    for (int i=0; i<numrates; i++){
sofia@0:      oindex = ii*ratenumf[i]*gpuarrays.osize[i];
sofia@0:      cudaDeviceSynchronize(); 
sofia@0:      memcpy(&h_out[i][oindex], gpuarrays.h_out[i], gpuarrays.osize[i]* ratenumf[i] *sizeof(float) );
sofia@0:    }
sofia@0:  }
sofia@0: 
sofia@0:  if (args.tim){
sofia@0:    cudaDeviceSynchronize();
sofia@0:    gettimeofday(&t_end, NULL); 
sofia@0:    t_gpu = (double) (t_end.tv_sec + (t_end.tv_usec / 1000000.0)  - t_start.tv_sec - (t_start.tv_usec/ 1000000.0)) * 1000.0; 
sofia@0:    printf("Finished GPU FIR\nProcessing gpu took: %f ms\n", t_gpu/(double)N);   
sofia@0:  }
sofia@0: 
sofia@0:  // check results
sofia@0:  check_results(h_reference, h_out, &gpuarrays, &gparams, N);
sofia@0:   
sofia@0:  //close GPU
sofia@0:  cudaMultiFilterFirClose(&gparams, &gpuarrays);
sofia@0:  
sofia@0:  // cleanup memory
sofia@0:  for (int i = 0; i < numrates; ++i){
sofia@0:    free(h_out[i]);
sofia@0:    free(h_in[i]);
sofia@0:    free(h_reference[i]);
sofia@0:  }
sofia@0:  
sofia@0:  free(farr.bk);
sofia@0:  
sofia@0: }