// SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
// SPDX-License-Identifier: BSD-3-Clause
#include "vtkPointGaussianVS.h"

const char *vtkPointGaussianVS =
"//VTK::System::Dec\n"
"\n"
"// SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen\n"
"// SPDX-License-Identifier: BSD-3-Clause\n"
"// this shader implements imposters in OpenGL for Spheres\n"
"\n"
"uniform float scaleFactor;\n"
"uniform int cameraParallel;\n"
"\n"
"// low-pass filtering matrix\n"
"// stored as a vector since it's a 2x2 symmetric matrix\n"
"uniform vec3 lowpassMatrix;\n"
"\n"
"in vec4 vertexMC;\n"
"\n"
"//VTK::Covariance::Dec\n"
"\n"
"out mat2 transformVCVSOutput;\n"
"\n"
"// optional normal declaration\n"
"//VTK::Normal::Dec\n"
"\n"
"// Texture coordinates\n"
"//VTK::TCoord::Dec\n"
"\n"
"// material property values\n"
"//VTK::Color::Dec\n"
"\n"
"// clipping plane vars\n"
"//VTK::Clip::Dec\n"
"\n"
"// camera and actor matrix values\n"
"//VTK::Camera::Dec\n"
"\n"
"// picking support\n"
"//VTK::Picking::Dec\n"
"\n"
"mat3 quaternionToMatrix(vec4 quat)\n"
"{\n"
"  // Normalize quaternion (required?)\n"
"  vec4 q = quat / length(quat);\n"
"  float r = q.x;\n"
"  float x = q.y;\n"
"  float y = q.z;\n"
"  float z = q.w;\n"
"\n"
"  return mat3(1. - 2. * (y * y + z * z), 2. * (x * y + r * z), 2. * (x * z - r * y),\n"
"    2. * (x * y - r * z), 1. - 2. * (x * x + z * z), 2. * (y * z + r * x),\n"
"    2. * (x * z + r * y), 2. * (y * z - r * x), 1. - 2. * (x * x + y * y));\n"
"}\n"
"\n"
"mat3 computeCov3D(vec3 scale, vec4 rotation)\n"
"{\n"
"  // Create scaling matrix\n"
"  mat3 S = mat3(\n"
"    scaleFactor * scale.x, 0., 0.,\n"
"    0., scaleFactor * scale.y, 0.,\n"
"    0., 0., scaleFactor * scale.z);\n"
"\n"
"  // Compute rotation matrix from quaternion\n"
"  mat3 R = quaternionToMatrix(rotation);\n"
"\n"
"  mat3 M = R * S;\n"
"\n"
"  // Compute 3D world covariance matrix\n"
"  return M * transpose(M);\n"
"}\n"
"\n"
"mat3 getProjectionMatrix(vec3 positionVC)\n"
"{\n"
"  // orthographic projection\n"
"  if (cameraParallel != 0)\n"
"  {\n"
"    return mat3(\n"
"      VCDCMatrix[0][0], 0., 0.,\n"
"      0., VCDCMatrix[1][1], 0.,\n"
"      0., 0., 0.);\n"
"  }\n"
"\n"
"  // affine approximation of the perspective projection\n"
"  // Refer to eq 29 and 31 in \"EWA Splatting\" (Zwicker et al., 2002).\n"
"  return mat3(\n"
"    VCDCMatrix[0][0] / positionVC.z, 0., 0.,\n"
"    0., VCDCMatrix[1][1] / positionVC.z, 0.,\n"
"    -(VCDCMatrix[0][0] * positionVC.x) / (positionVC.z * positionVC.z), -(VCDCMatrix[1][1] * positionVC.y) / (positionVC.z * positionVC.z), 0.);\n"
"}\n"
"\n"
"void main()\n"
"{\n"
"  //VTK::Color::Impl\n"
"\n"
"  //VTK::Normal::Impl\n"
"\n"
"  //VTK::TCoord::Impl\n"
"\n"
"  //VTK::Clip::Impl\n"
"\n"
"  vec4 posVC = MCVCMatrix * vertexMC;\n"
"\n"
"  mat3 J = getProjectionMatrix(posVC.xyz);\n"
"\n"
"  mat3 W = mat3(MCVCMatrix);\n"
"\n"
"  mat3 T = J * W;\n"
"\n"
"  //VTK::Covariance::Impl\n"
"\n"
"  // discard third column/row and store the matrix in a vector\n"
"  mat2 cov2d = mat2(cov);\n"
"\n"
"  // It is possible to apply a convolution here by adding a custom 2D covariance matrix\n"
"  // For example, it can be useful for making sure the point is at least 1 pixel wide\n"
"  cov2d += mat2(lowpassMatrix.x, lowpassMatrix.y, lowpassMatrix.y, lowpassMatrix.z);\n"
"\n"
"  // eigen vectors of the covariance matrix\n"
"  if (abs(cov2d[0][1]) > 1e-6)\n"
"  {\n"
"    // compute basis transformation based on eigen decomposition of the 2D covariance\n"
"    // half of the covariance matrix trace\n"
"    float halfTrace = 0.5 * (cov2d[0][0] + cov2d[1][1]);\n"
"\n"
"    // eigen values of the covariance matrix\n"
"    float term = sqrt(halfTrace * halfTrace - determinant(cov2d));\n"
"    float eigenValue1 = halfTrace + term;\n"
"    float eigenValue2 = halfTrace - term;\n"
"\n"
"    vec2 eigenVector1 = normalize(vec2(eigenValue1 - cov2d[1][1], cov2d[0][1]));\n"
"    vec2 eigenVector2 = vec2(eigenVector1.y, -eigenVector1.x);\n"
"\n"
"    transformVCVSOutput = mat2(sqrt(eigenValue1) * eigenVector1, sqrt(eigenValue2) * eigenVector2);\n"
"  }\n"
"  else\n"
"  {\n"
"    // the covariance is a diagonal matrix, so axis aligned\n"
"    // we just need to extract the scaling matrix (the square root of the covariance in this case)\n"
"    transformVCVSOutput = mat2(sqrt(cov2d[0][0]), 0.0, 0.0, sqrt(cov2d[1][1]));\n"
"  }\n"
"\n"
"  // filter out points outside of the frustum inflated by 30% because the approximation may stretch them too much\n"
"  // see https://github.com/graphdeco-inria/gaussian-splatting/issues/191\n"
"  // this is a hand-tuned parameter but it should work on most cases (if not, it may be useful to expose the value as uniform)\n"
"  vec4 posDC = VCDCMatrix * posVC;\n"
"\n"
"  if (posDC.x < -1.3*posDC.w || posDC.x > 1.3*posDC.w || posDC.y < -1.3*posDC.w || posDC.y > 1.3*posDC.w)\n"
"  {\n"
"    transformVCVSOutput = mat2(0);\n"
"  }\n"
"\n"
"  gl_Position = posDC;\n"
"\n"
"  //VTK::Picking::Impl\n"
"}\n"
"";