0

I am trying to make a zoom_image function which zooms gray image using discrete Fourier transform,. the code i am include work if the image size is less than or equal 4*4 but if size increase. it gives 'double free or corruption (out) Aborted (core dumped)' error

I have tried fft_d and ifft_2d function of my code it works if input size is small if input size is big it gives the above error

#include<math.h>
#include<stdio.h>
#include<stdlib.h>
struct complex{
      float real;
      float im;
};
typedef struct complex complex;
complex w(int i, int n) {
      complex result;
      result.real = cos(2*M_PI*i/n);
      result.im = -sin(2*M_PI*i/n);
      return result;
}
complex wp(int i, int n) {
      complex result;
      result.real = cos(2*M_PI*i/n);
      result.im = sin(2*M_PI*i/n);
      return result;
}
complex mul(complex a, complex b) {
      complex result;
      result.real = a.real*b.real - a.im*b.im;
      result.im = a.real*b.im + b.real*a.im;
      return result;
}
complex divi(complex a, complex b) {
      complex result;
      result.real = (a.real*b.real + a.im*b.im)/(b.real*b.real + b.im*b.im);
      result.im = (-a.real*b.im + b.real*a.im)/(b.real*b.real + b.im*b.im);
      return result;
}
complex add(complex a, complex b) {
      complex result;
      result.real = a.real+b.real ;
      result.im = a.im+b.im;
      return result;
}
complex sub(complex a, complex b) {
      complex result;
      result.real = a.real - b.real;
      result.im = a.im - b.im;
      return result;
}
void printComplex(complex var) {
      printf("%f i%f\n",var.real,var.im);
}
void printComplexS(complex var) {
      printf("%f i%f",var.real,var.im);
}

#include<stdio.h>
#include<stdlib.h>
#include"complex.h"
static int gb=0;
complex* _fft(complex *arr, int size, int step, int index) {
    if(size == 2) {
        complex *result;
        result = (complex*)malloc(size*sizeof(complex));
        result[0] = add(arr[index], arr[index+step]);
        result[1] = sub(arr[index], arr[index+step]);
        return result;
    }
    else {
        int i;
        complex *even, *odd, *result, *mull;
        even = _fft(arr, size/2, step*2, index);
        odd = _fft(arr, size/2, step*2, index+step);
        result = (complex*)malloc(size*sizeof(complex));
        mull = (complex*)malloc(size/2*sizeof(complex));
        for(i=0;i<size/2;i++) {
            mull[i] = mul(odd[i], w(i,size));
        }
        for(i=0;i<size/2;i++)
            result[i] = add(even[i], mull[i]);
        for(;i<size;i++)
            result[i] = sub(even[i - size/2], mull[i - size/2]);
        free(even);
        free(odd);
        free(mull);
        return result;
    }
}
complex* fft(complex *arr, int size) {
      return (complex*)_fft(arr, size, 1, 0);
}
complex* _ifft(complex *arr, int size, int step, int index) {
      if(size == 2) {
            complex *result;
            result = (complex*)malloc(size*sizeof(complex));
            result[0] = add(arr[index], arr[index+step]);
            result[1] = sub(arr[index], arr[index+step]);
            return result;
      }
      else {
            int i;
            complex *even, *odd, *result, *mull;
            even = _ifft(arr, size/2, step*2, index);
            odd = _ifft(arr, size/2, step*2, index+step);
            result = (complex*)malloc(size*sizeof(complex));
            mull = (complex*)malloc(size/2*sizeof(complex));
            for(i=0;i<size/2;i++)
                  mull[i] = mul(odd[i], wp(i,size));
            for(i=0;i<size/2;i++)
                  result[i] = add(even[i], mull[i]);
            for(;i<size;i++)
                  result[i] = sub(even[i - size/2], mull[i - size/2]);
            free(even);
            free(odd);
            free(mull);
            return result;
      }
}
complex* ifft(complex *arr, int size) {
    complex *re = _ifft(arr, size, 1, 0);
    for(int i=0;i<size;i++){
        re[i].real=re[i].real/size;
        re[i].im=re[i].im/size;
    }
    return re;
}
complex** transpose(complex **src, int n, int m) {
    complex **result, ***re;
    result=(complex**)malloc(sizeof(complex*)*m);
    for(int i=0;i<m;i++) {
        result[i]=(complex*)malloc(sizeof(complex)*n);
        for(int j=0;j<n;j++)
            result[i][j]=src[j][i];
    }
    re=(complex***)&result;
    return (complex**)*re;
}
complex** fft_2d(complex** arr, int n, int m) {
    complex **arrR, ***r, *temp;
    int i, j;
    arrR = (complex**)malloc(sizeof(complex*));
    for(i=0;i<n;i++) {
        arrR[i]=(complex*)fft(arr[i],m);
        printf("%d ",i);
    }
    arrR=(complex**)transpose(arrR,n,m);
    malloc(0);
    for(i=0;i<m;i++)
        arrR[i]=(complex*)fft(arrR[i],n);
    arrR=transpose(arrR,n,m);
    r=(complex***)&arrR;
    return (complex**)*r;
}
complex** ifft_2d(complex** arr, int n, int m) {
    complex **arrR, ***r, *temp;
    int i, j;
    arrR = (complex**)malloc(sizeof(complex*));
    for(i=0;i<n;i++)
        arrR[i]=(complex*)ifft(arr[i],m);
    arrR=(complex**)transpose(arrR,n,m);
    malloc(0);
    for(i=0;i<m;i++)
        arrR[i]=(complex*)ifft(arrR[i],n);
    arrR=transpose(arrR,n,m);
    r=(complex***)&arrR;
    return (complex**)*r;
}

unsigned int** zoom_img(unsigned int **img, int col, int row, int r_col, int r_row) {
    int i, j;
    complex **mat=(complex**)malloc(sizeof(complex*)*col), **re;
    complex **new_re=(complex**)malloc(sizeof(complex*)*(col+2*r_col)), **u;
    unsigned int **result=(unsigned int**)malloc(sizeof(unsigned int*)*(col + 2*r_col)), ***z;
    for(i=0;i<col;i++)
        mat[i]=(complex*)malloc(sizeof(complex)*row);
    for (i=0;i<col;i++) {
        for (j=0;j<row;j++) {
            mat[i][j].real = (float)pow(-1, i+j)*(float)img[i][j];
            mat[i][j].im=0;
        }
    }
    re = (complex**)fft_2d(mat, col, row);
    for(i=0;i<(col+2*r_col);i++)
        new_re[i]=(complex*)malloc(sizeof(complex)*(row+2*r_row));
    for(i=0;i<(col+2*r_col);i++) {
        for(j=0;j<(row+2*r_row);j++) {
            if(i<r_col || i>r_col+col-1 || j<r_row || j>r_row+row-1) {
                new_re[i][j].real = 0;
                new_re[i][j].im = 0;
            }
            else
                new_re[i][j]=re[i-r_col][j-r_row];
        }
    }
    u = (complex**)ifft_2d(new_re, col+2*r_col, row + 2*r_row);
    for(i=0;i<(col+2*r_col);i++) {
        result[i]=(unsigned int*)malloc(sizeof(unsigned int)*(row+2*r_row));
        for(j=0;j<(row+2*r_row);j++) {
            result[i][j] = (unsigned int)u[i][j].real;
        }
    }
    z=&result;
    return *z;
}

int main() {
    unsigned int i, j, **arr=(unsigned int**)malloc(sizeof(unsigned int*)*2), **result;
    for(i=0;i<2;i++) {
        arr[i]=(unsigned int*)malloc(sizeof(unsigned int)*2);
        for(j=0;j<2;j++) {
             arr[i][j] = i+j +1;
        }
    }
    result = zoom_img(arr,2,2,2,2);
    return 0;
}/*
int main() {
    complex **arr, **result, **re;
    arr=(complex**)malloc(sizeof(complex*)*4);
    for (int i = 0; i < 4; i++) {
        arr[i]=(complex*)malloc(sizeof(complex)*4);
        for (int j = 0; j < 4; j++) {
            arr[i][j].real = i*j+1.2;
            arr[i][j].im=0;
        }
    }
    result = (complex**)fft_2d(arr,4,4);
    //malloc(0);
    re = (complex**)ifft_2d(result,4,4);
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 4; j++) {
            printf("(%2f) ", arr[i][j].real);
            printComplexS(result[i][j]);
            printf(" (%2f) ",re[i][j].real);
        }
        printf("\n");
    }
      return 0;
}*/
Kishan Dhrangadhariya
  • 3
  • 2
  • Is `size` guaranteed to be a power of two? Otherwise _fft had issues with boundary conditions, where it ends up being called with size=1, which leads to calls with size=1/2, I.e. size=0, and malloc(0) – Michael Veksler Mar 30 '19 at 18:24
  • After adding print statements, it is evident that size is not a power of two, which leads to the crash – Michael Veksler Mar 30 '19 at 18:51
  • size would never reduce to 1 because it when size is 2 it will directly go to if statement – Kishan Dhrangadhariya Mar 31 '19 at 08:54
  • and yes size is guaranteed to be in power of 2 – Kishan Dhrangadhariya Mar 31 '19 at 08:54
  • In that case I suggest adding an assert that validates that the value is in fact a power of two, since with you `main` a value of 6 is passed – Michael Veksler Mar 31 '19 at 09:43
  • I tried using fft_2d alone on my system with size of 32 it gives double-free or corruption (out) error – Kishan Dhrangadhariya Mar 31 '19 at 10:27
  • so please post that code instead, since `zoom_img` makes calls which are guaranteed not to be powers of two: `u = (complex**)ifft_2d(new_re, col+2*r_col, row + 2*r_row);` – Michael Veksler Mar 31 '19 at 10:49
  • Side note: FFT is usually done with a couple of arrays allocated at the beginning and reused all the time. The way you use malloc (with so many calls) will make your code run much slower than it could have. I have implemented FFT too many years ago, so I don't remember the details, but you can look it up – Michael Veksler Mar 31 '19 at 10:55

1 Answers1

0

The code has several issues. 

complex** ifft_2d(complex** arr, int n, int m) {
    complex **arrR, ***r;
    int i;

And then you are allocating too little: 

    arrR = (complex**)malloc(sizeof(complex*));

The above line should be malloc(sizeof(complex*) * n). Without the *n the next line causes undefined behavior (buffer overflow):

    for(i=0;i<n;i++)
        arrR[i]=(complex*)ifft(arr[i],m);

Then there is this strange line, why is it for?

    malloc(0);

Then the ending of the function is weird (Why not simply return arrR?):

    r=(complex***)&arrR;
    return (complex**)*r;

Next, the recursive code (_ifft and _fft) assumes that size is a power of two, otherwise it gets into infinite recursion. You should validate the input. At least assert that:

assert(n >= 2);
assert(((n-1) & n ) == 0); // for positive n, assert that it is a power of 2

With this line you can see that zoom_img violates this assumption in the line:

u = (complex**)ifft_2d(new_re, col+2*r_col, row + 2*r_row);

Note that in the uncommented main, col==2, row==2, r_col==2, r_row==2, which ends up with size == 6, which is not a power of 2.

A smaller issue is performance. The code overuses malloc instead of reusing the same block of memory, and peek into different areas of it. This is what the classical FFT does. Classical FFT also does not use recursion like that, but uses iterations instead.

Michael Veksler
  • 8,217
  • 1
  • 20
  • 33