Implementation of dct forward and inverse

dct.c

#include <math.h>

void dct(float *dct_matrix, float *matrix, int N){
  int i,j,u,v;
  float Cu, Cv;
  memset(dct_matrix, 0, N*N*sizeof(float));
  float accum = 0;
  for (u = 0; u < N; ++u){
    for (v = 0; v < N; ++v){
      if (u == 0){
        Cu = 1.0/sqrt(2.0);
      }else{
        Cu = 1.0;
      }
      if (v == 0){
        Cv = 1.0 /sqrt(2.0);
      }else{
        Cv = 1.0;
      }
      accum = 0.0f;
      for (i = 0; i < N; i ++){
        for (j = 0; j < N; j++){
          accum += matrix[i*8+j] * cos((2*j+1)*u*M_PI/16.0) * cos((2*i+1)*v*M_PI/16.0);
        }
      }
      dct_matrix[u*8+v] = 0.25*Cu*Cv*accum;
    }
  }
}

void idct(float *matrix, float *dct_matrix, int N){
  int i,j,u,v;
  float Cu, Cv;
  float accum = 0;
  for (i = 0; i < N; ++i){
    for (j = 0; j < N; ++j){
      accum = 0.0f;
      for (u = 0; u < N; u ++){
        for (v = 0; v < N; v++){
          if (u == 0){
            Cu = 1.0/sqrt(2.0);
          }else{
            Cu = 1.0;
          }
          if (v == 0){
            Cv = 1.0 /sqrt(2.0);
          }else{
            Cv = 1.0;
          }
          accum += dct_matrix[u*8+v] * Cu * Cv * cos((2*j+1)*u*M_PI/16.0) * cos((2*i+1)*v*M_PI/16.0);
        }
      }
      matrix[i*8+j] = 0.25*accum;
    }
  }
}