

int N,M;
int *OFFSET;
void setOffset(){
  int W,H,X[4];
  int i,j,cnt=0;
  
  W=(IMG->width-21)/236+1;
  H=(IMG->height-21)/236+1;
  N=W*H;
  M=(N-1)/2+1;
  OFFSET=malloc(sizeof(int)*N);
  X[0]=236;
  X[1]=236*IMG->width;
  X[2]=IMG->width-256;
  X[3]=IMG->width*(IMG->height-256);
  for(j=0;j<H-1;j++){
    for(i=0;i<W-1;i++){
      OFFSET[cnt]=X[0]*i+X[1]*j;
      cnt++;
    }
  }
  for(j=0;j<H-1;j++){
    OFFSET[cnt]=X[2]+X[1]*j;
    cnt++;
  }
  for(i=0;i<W-1;i++){
    OFFSET[cnt]=X[0]*i+X[3];
    cnt++;
  }
  OFFSET[cnt]=X[2]+X[3];
}

void cropBlock(
int n,
double *io){
  int i,j;
  
  for(j=0;j<256;j++){
    for(i=0;i<256;i++){
      fft_in[i+256*j][0]=io[OFFSET[n]+i+j*IMG->width];
      fft_in[i+256*j][1]=io[OFFSET[N-n-1]+i+j*IMG->width];
    }
  }
}

void embedBlock(
int n,
double *io){
  int i,j;
  
  for(j=10;j<246;j++){
    for(i=10;i<246;i++){
      io[OFFSET[n]+i+j*IMG->width]=fft_out[i+256*j][0];
      io[OFFSET[N-n-1]+i+j*IMG->width]=fft_out[i+256*j][1];
    }
  }
}

void exeFilterTile(
double *io){
  int i,n;
  double *tmp;
  
  tmp=malloc(sizeof(double)*IMG->width*IMG->height);
  for(i=0;i<IMG->width*IMG->height;i++){
    tmp[i]=io[i];
  }
  for(n=0;n<M;n++){
    cropBlock(n,tmp);
    fftw_execute(p_f);
    for(i=0;i<65536;i++){
      fft_in[i][0]=fft_blur[i][0]*fft_out[i][0]-fft_blur[i][1]*fft_out[i][1];
      fft_in[i][1]=fft_blur[i][0]*fft_out[i][1]+fft_blur[i][1]*fft_out[i][0];
    }
    fftw_execute(p_i);
    embedBlock(n,io);
  }
  free(tmp);
}

int main(){
  clock_t start,stop;
  
  defBlur();
  readImage("input.bmp");
  splitRGB();
  start=clock();
  setFFT(256,256);
  transFilter(256,256);
  setOffset();
  exeFilterTile(R);
  exeFilterTile(G);
  exeFilterTile(B);
  freeFFT();
  free(OFFSET);
  stop=clock();
  stop-=start;
  mergeRGB();
  writeImage("output.bmp");
  printf("%d\n",stop);
  getch();
  return 0;
}
FIR フィルタリングを FFT で高速化してみよう
