#include "ccv.h"

#include "ccv_internal.h"

void ccv_decimal_slice(ccv_dense_matrix_t* a, ccv_dense_matrix_t** b, int type, float y, float x, int rows, int cols)

{

  ccv_declare_derived_signature(sig, a->sig != 0, ccv_sign_with_format(128, "ccv_decimal_slice(%a,%a,%d,%d)", y, x, rows, cols), a->sig, CCV_EOF_SIGN);

  type = (type == 0) ? CCV_GET_DATA_TYPE(a->type) | CCV_GET_CHANNEL(a->type) : 
CCV_GET_DATA_TYPE0
(type) | 0
CCV_GET_CHANNEL0
(a->type);

  ccv_dense_matrix_t* db = *b = ccv_dense_matrix_renew(*b, rows, cols, CCV_ALL_DATA_TYPE | CCV_GET_CHANNEL(a->type), type, sig);

  ccv_object_return_if_cached(, db);

  int i, j, ch = CCV_GET_CHANNEL(a->type);

  int ix = (int)x, iy = (int)y;

  float xd = x - ix, yd = y - iy;

  float w00 = (1 - xd) * (1 - yd);

  float w01 = xd * (1 - yd);

  float w10 = (1 - xd) * yd;

  float w11 = xd * yd;

  int dx = 0, dy = 0;

  int rows_1 = 0, cols_1 = 0;

  // it is going to be hard to deal with border efficiently, since this is used for tld, will ignore border for now

  if (!(iy >= 0 && iy + rows < a->rows && ix >= 0 && ix + cols < a->cols))

  {

    ccv_zero(db);

    if (iy < 0) { rows += iy; dy = -iy; iy = 0; }

    if (iy + rows >= a->rows)

    {

      rows = a->rows - iy - 1;

      if (iy + rows > a->rows)

        rows_1 = 1; // we need to do our best to padding the last row

    }

    if (ix < 0) { cols += ix; dx = -ix; ix = 0; }

    if (x + cols >= a->cols)

    {

      cols = a->cols - ix - 1;

      if (x + cols > a->cols)

        cols_1 = 1; // we need to do our best to padding the last col

    }

  }

  unsigned char* a_ptr = (unsigned char*)ccv_get_dense_matrix_cell(a, iy, ix, 0);

  unsigned char* b_ptr = (unsigned char*)ccv_get_dense_matrix_cell(db, dy, dx, 0);

#define for_block(_for_set, _for_get) \

  
for (i = 0; 1
i < rows; 
i++111
) \

  { \

    for (j = 0; j < cols * ch; 
j++30.3k
) \

    { \

      _for_set(b_ptr, j, (_for_get(a_ptr, j) * G00 + _for_get(a_ptr, j + ch) * G01 + _for_get(a_ptr + a->step, j) * G10 + _for_get(a_ptr + a->step, j + ch) * G11) / GALL); \

    } \

    if (cols_1) \

      _for_set(b_ptr, j, (_for_get(a_ptr, j) * (G00 + G01) + _for_get(a_ptr + a->step, j) * G10 + _for_get(a_ptr + a->step, j + ch) * G11) / GALL); \

    a_ptr += a->step; \

    b_ptr += db->step; \

  } \

  if (rows_1) \

  { \

    for (j = 0; j < cols * ch; j++) \

    { \

      _for_set(b_ptr, j, (_for_get(a_ptr, j) * (G00 + G10) + _for_get(a_ptr, j + ch) * (G01 + G11)) / GALL); \

    } \

    if (cols_1) \

      _for_set(b_ptr, j, _for_get(a_ptr, j)); \

  }

  /* unswitch in the manual way so that we can use integer interpolation */

  if ((a->type & CCV_8U) || 
(a->type & CCV_32S)0
 || 
(a->type & CCV_64S)0
)

  {

    const int W_BITS14 = 14;

    int iw00 = (int)(w00 * (1 << W_BITS14) + 0.5);

    int iw01 = (int)(w01 * (1 << W_BITS14) + 0.5);

    int iw10 = (int)(w10 * (1 << W_BITS14) + 0.5);

    int iw11 = (1 << W_BITS14) - iw00 - iw01 - iw10;

#define G00 (iw00)

#define G01 (iw01)

#define G10 (iw10)

#define G11 (iw11)

#define GCOM (1 << (W_BITS14 - 1))

#define GALL (1 << (W_BITS14))

    ccv_matrix_setter(db->type, ccv_matrix_getter_integer_only, a->type, for_block);

#undef G00

#undef G01

#undef G10

#undef G11

#undef GCOM

#undef GALL

  } else {

#define G00 (w00)

#define G01 (w01)

#define G10 (w10)

#define G11 (w11)

#define GCOM (0)

#define GALL (1)

    ccv_matrix_setter(db->type, ccv_matrix_getter_float_only, a->type, for_block);

#undef G00

#undef G01

#undef G10

#undef G11

#undef GCOM

#undef GALL

  }

#undef for_block

}

ccv_decimal_point_t ccv_perspective_transform_apply(ccv_decimal_point_t point, ccv_size_t size, float m00, float m01, float m02, float m10, float m11, float m12, float m20, float m21, float m22)

{

  m00 *= 1.0 / ccv_max(size.width, size.height);

  m01 *= 1.0 / ccv_max(size.width, size.height);

  m02 *= 1.0 / ccv_max(size.width, size.height);

  m10 *= 1.0 / ccv_max(size.width, size.height);

  m11 *= 1.0 / ccv_max(size.width, size.height);

  m12 *= 1.0 / ccv_max(size.width, size.height);

  m20 *= 1.0 / (ccv_max(size.width, size.height) * ccv_max(size.width, size.height));

  m21 *= 1.0 / (ccv_max(size.width, size.height) * ccv_max(size.width, size.height));

  m22 *= 1.0 / ccv_max(size.width, size.height);

  point.x -= size.width * 0.5;

  point.y -= size.height * 0.5;

  float wz = 1.0 / (point.x * m20 + point.y * m21 + m22);

  float wx = size.width * 0.5 + (point.x * m00 + point.y * m01 + m02) * wz;

  float wy = size.height * 0.5 + (point.x * m10 + point.y * m11 + m12) * wz;

  return ccv_decimal_point(wx, wy);

}

// this method is a merely baseline implementation and has no optimization effort ever put into it, if at all

void ccv_perspective_transform(ccv_dense_matrix_t* a, ccv_dense_matrix_t** b, int type, float m00, float m01, float m02, float m10, float m11, float m12, float m20, float m21, float m22)

{

  ccv_declare_derived_signature(sig, a->sig != 0, ccv_sign_with_format(128, "ccv_perspective_transform(%a,%a,%a,%a,%a,%a,%a,%a,%a)", m00, m01, m02, m10, m11, m12, m20, m21, m22), a->sig, CCV_EOF_SIGN);

  type = (type == 0) ? CCV_GET_DATA_TYPE(a->type) | CCV_GET_CHANNEL(a->type) : 
CCV_GET_DATA_TYPE0
(type) | 0
CCV_GET_CHANNEL0
(a->type);

  ccv_dense_matrix_t* db = *b = ccv_dense_matrix_renew(*b, a->rows, a->cols, CCV_ALL_DATA_TYPE | CCV_GET_CHANNEL(a->type), type, sig);

  ccv_object_return_if_cached(, db);

  // with default of bilinear interpolation

  int i, j, k, ch = CCV_GET_CHANNEL(a->type);

  unsigned char* a_ptr = a->data.u8;

  unsigned char* b_ptr = db->data.u8;

  // assume field of view is 60, modify the matrix value to reflect that

  // (basically, apply x / ccv_max(a->rows, a->cols), y / ccv_max(a->rows, a->cols) before hand

  m00 *= 1.0 / ccv_max(a->rows, a->cols);

  m01 *= 1.0 / ccv_max(a->rows, a->cols);

  m02 *= 1.0 / ccv_max(a->rows, a->cols);

  m10 *= 1.0 / ccv_max(a->rows, a->cols);

  m11 *= 1.0 / ccv_max(a->rows, a->cols);

  m12 *= 1.0 / ccv_max(a->rows, a->cols);

  m20 *= 1.0 / (ccv_max(a->rows, a->cols) * ccv_max(a->rows, a->cols));

  m21 *= 1.0 / (ccv_max(a->rows, a->cols) * ccv_max(a->rows, a->cols));

  m22 *= 1.0 / ccv_max(a->rows, a->cols);

#define for_block(_for_set, _for_get) \

  
for (i = 0; 1
i < db->rows; 
i++500
) \

  { \

    float cy = i - db->rows * 0.5; \

    float crx = cy * m01 + m02; \

    float cry = cy * m11 + m12; \

    float crz = cy * m21 + m22; \

    for (j = 0; j < db->cols; 
j++250k
) \

    { \

      float cx = j - db->cols * 0.5; \

      float wz = 1.0 / (cx * m20 + crz); \

      float wx = a->cols * 0.5 + (cx * m00 + crx) * wz; \

      float wy = a->rows * 0.5 + (cx * m10 + cry) * wz; \

      int iwx = (int)wx; \

      int iwy = (int)wy; \

      wx = wx - iwx; \

      wy = wy - iwy; \

      int iwx1 = ccv_min(iwx + 1, a->cols - 1); \

      int iwy1 = ccv_min(iwy + 1, a->rows - 1); \

      if (iwx >= 0 && iwx < a->cols && 
iwy >= 0222k
 && 
iwy < a->rows208k
) \

        
for (k = 0; 195k
k < ch; 
k++586k
) \

          _for_set(b_ptr, j * ch + k, _for_get(a_ptr + iwy * a->step, iwx * ch + k) * (1 - wx) * (1 - wy) + \

                        _for_get(a_ptr + iwy * a->step, iwx1 * ch + k) * wx * (1 - wy) + \

                        _for_get(a_ptr + iwy1 * a->step, iwx * ch + k) * (1 - wx) * wy + \

                        _for_get(a_ptr + iwy1 * a->step, iwx1 * ch + k) * wx * wy); \

      else \

        
for (k = 0; 54.5k
k < ch218k
; 
k++163k
) \

          _for_set(b_ptr, j * ch + k, 0); \

    } \

    b_ptr += db->step; \

  }

  ccv_matrix_setter(db->type, ccv_matrix_getter, a->type, for_block);

#undef for_block

}