#include "projection_constraint.h"

template <
  typename DerivedUV,
  typename DerivedM,
  typename DerivedVP,
  typename DerivedA,
  typename DerivedB>
void igl::projection_constraint(
  const Eigen::MatrixBase<DerivedUV> & UV,
  const Eigen::MatrixBase<DerivedM> & _M,
  const Eigen::MatrixBase<DerivedVP> & VP,
  Eigen::PlainObjectBase<DerivedA> & A,
  Eigen::PlainObjectBase<DerivedB> & B)
{
  typedef typename DerivedA::Scalar Scalar;
  const Scalar u = UV(0);
  const Scalar v = UV(1);
  const Scalar cu = VP(0);
  const Scalar cv = VP(1);
  const Scalar w = VP(2);
  const Scalar h = VP(3);
  // u = cu + w*(0.5 + 0.5*((M.row(0)*X) / (M.row(3)*X) ))
  // u-cu = w*(0.5 + 0.5*((M.row(0)*X) / (M.row(3)*X) ))
  // (u-cu)/w = 0.5 + 0.5*((M.row(0)*X) / (M.row(3)*X) )
  // (u-cu)/w - 0.5 = 0.5*((M.row(0)*X) / (M.row(3)*X) )
  // 2.*(u-cu)/w - 1 = ((M.row(0)*X) / (M.row(3)*X) )
  // (2.*(u - cu)/w - 1) * M.row(3)*X = M.row(0)*X
  // (2.*(u - cu)/w - 1) * M.row(3)*X - M.row(0)*X = 0
  // (2.*(u - cu)/w - 1) * (M.block(3,0,1,3)*x + M(3,3)) - M.block(0,0,1,3)*x - M(0,3) = 0
  // (2.*(u - cu)/w - 1) * (M.block(3,0,1,3)*x + M(3,3)) - M.block(0,0,1,3)*x = M(0,3)
  // ((2.*(u - cu)/w - 1) * M.block(3,0,1,3) - M.block(0,0,1,3))*x = M(0,3) - (2.*(u - cu)/w - 1)*M(3,3)
  Eigen::Matrix<Scalar,4,4> M = _M.template cast<Scalar>();
  A.resize(2,3);
  A<<
   ((2.*(u - cu)/w - 1.) * M.block(3,0,1,3) - M.block(0,0,1,3)),
   ((2.*(v - cv)/h - 1.) * M.block(3,0,1,3) - M.block(1,0,1,3));
  B.resize(2,1);
  B<<
   M(0,3) - (2.*(u - cu)/w - 1.)*M(3,3),
   M(1,3) - (2.*(v - cv)/h - 1.)*M(3,3);
}

#ifdef IGL_STATIC_LIBRARY
// Explicit template instantiation
// generated by autoexplicit.sh
template void igl::projection_constraint<Eigen::Matrix<double, -1, 1, 0, 3, 1>, Eigen::Matrix<float, 4, 4, 0, 4, 4>, Eigen::Matrix<float, 4, 1, 0, 4, 1>, Eigen::Matrix<double, 2, 3, 0, 2, 3>, Eigen::Matrix<double, 2, 1, 0, 2, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, 3, 1> > const&, Eigen::MatrixBase<Eigen::Matrix<float, 4, 4, 0, 4, 4> > const&, Eigen::MatrixBase<Eigen::Matrix<float, 4, 1, 0, 4, 1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 3, 0, 2, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> >&);
// generated by autoexplicit.sh
template void igl::projection_constraint<Eigen::Matrix<double, 2, 1, 0, 2, 1>, Eigen::Matrix<float, 4, 4, 0, 4, 4>, Eigen::Matrix<float, 4, 1, 0, 4, 1>, Eigen::Matrix<double, 2, 3, 0, 2, 3>, Eigen::Matrix<double, 2, 1, 0, 2, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> > const&, Eigen::MatrixBase<Eigen::Matrix<float, 4, 4, 0, 4, 4> > const&, Eigen::MatrixBase<Eigen::Matrix<float, 4, 1, 0, 4, 1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 3, 0, 2, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> >&);
#endif