/*=========================================================================

  Program:   Visualization Toolkit
  Module:    BoxClipTriangulate.cxx

  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
  All rights reserved.
  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
     PURPOSE.  See the above copyright notice for more information.

=========================================================================*/

/*----------------------------------------------------------------------------
 Copyright (c) Sandia Corporation
 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
----------------------------------------------------------------------------*/

// This test makes sure that vtkBoxClipDataSet correctly triangulates all cell
// types.

#include "vtkBoxClipDataSet.h"
#include "vtkCellArray.h"
#include "vtkDataSetSurfaceFilter.h"
#include "vtkDoubleArray.h"
#include "vtkMath.h"
#include "vtkPoints.h"
#include "vtkPolyData.h"
#include "vtkTriangle.h"
#include "vtkUnstructuredGrid.h"

#include "vtkActor.h"
#include "vtkExtractEdges.h"
#include "vtkPolyDataMapper.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkRenderer.h"

#include "vtkSmartPointer.h"
#define VTK_CREATE(type, var) vtkSmartPointer<type> var = vtkSmartPointer<type>::New()

#include <ctime>

#include <vector>

const int NumPoints = 13;
const double PointData[NumPoints * 3] = {
  0.0, 0.0, 0.0, //
  0.0, 1.0, 0.0, //
  1.0, 0.0, 0.0, //
  1.0, 1.0, 0.0, //
  2.0, 0.0, 0.0, //
  2.0, 1.0, 0.0, //

  0.0, 0.0, 1.0, //
  0.0, 1.0, 1.0, //
  1.0, 0.0, 1.0, //
  1.0, 1.0, 1.0, //
  2.0, 0.0, 1.0, //
  2.0, 1.0, 1.0, //
  2.0, 0.5, 1.0  //
};

const vtkIdType NumTriStripCells = 1;
const vtkIdType TriStripCells[] = {
  6, 1, 0, 3, 2, 5, 4 //
};

const vtkIdType NumQuadCells = 2;
const vtkIdType QuadCells[] = {
  4, 0, 2, 3, 1, //
  4, 2, 4, 5, 3  //
};

const vtkIdType NumPixelCells = 2;
const vtkIdType PixelCells[] = {
  4, 0, 2, 1, 3, //
  4, 2, 4, 3, 5  //
};

const vtkIdType NumPolyCells = 3;
const vtkIdType PolyCells[] = {
  4, 0, 2, 3, 1,    //
  3, 2, 4, 5,       //
  5, 6, 8, 12, 9, 7 //
};

const vtkIdType NumHexCells = 2;
const vtkIdType HexCells[] = {
  8, 6, 8, 2, 0, 7, 9, 3, 1,  //
  8, 4, 2, 8, 10, 5, 3, 9, 11 //
};
const vtkIdType NumExpectedHexSurfacePolys = 20;

const vtkIdType NumVoxelCells = 2;
const vtkIdType VoxelCells[] = {
  8, 0, 2, 1, 3, 6, 8, 7, 9,  //
  8, 10, 8, 11, 9, 4, 2, 5, 3 //
};
const vtkIdType NumExpectedVoxelSurfacePolys = 20;

const vtkIdType NumWedgeCells = 4;
const vtkIdType WedgeCells[] = {
  6, 0, 1, 2, 6, 7, 8,  //
  6, 7, 8, 9, 1, 2, 3,  //
  6, 8, 11, 9, 2, 5, 3, //
  6, 2, 5, 4, 8, 11, 10 //
};
const vtkIdType NumExpectedWedgeSurfacePolys = 20;

const vtkIdType NumPyramidCells = 2;
const vtkIdType PyramidCells[] = {
  5, 8, 9, 3, 2, 0, //
  5, 2, 3, 9, 8, 12 //
};
const vtkIdType NumExpectedPyramidSurfacePolys = 8;

class BoxClipTriangulateFailed
{
};

//-----------------------------------------------------------------------------

static void CheckWinding(vtkBoxClipDataSet* alg)
{
  alg->Update();
  vtkUnstructuredGrid* data = alg->GetOutput();

  vtkPoints* points = data->GetPoints();

  vtkCellArray* cells = data->GetCells();
  cells->InitTraversal();

  vtkIdType npts;
  const vtkIdType* pts;
  while (cells->GetNextCell(npts, pts))
  {
    if (npts != 4)
    {
      std::cout << "Weird.  I got something that is not a tetrahedra." << std::endl;
      continue;
    }

    double p0[3], p1[3], p2[3], p3[3];
    points->GetPoint(pts[0], p0);
    points->GetPoint(pts[1], p1);
    points->GetPoint(pts[2], p2);
    points->GetPoint(pts[3], p3);

    // If the winding is correct, the normal to triangle p0,p1,p2 should point
    // towards p3.
    double v0[3], v1[3];
    v0[0] = p1[0] - p0[0];
    v0[1] = p1[1] - p0[1];
    v0[2] = p1[2] - p0[2];
    v1[0] = p2[0] - p0[0];
    v1[1] = p2[1] - p0[1];
    v1[2] = p2[2] - p0[2];

    double n[3];
    vtkMath::Cross(v0, v1, n);

    double d[3];
    d[0] = p3[0] - p0[0];
    d[1] = p3[1] - p0[1];
    d[2] = p3[2] - p0[2];

    if (vtkMath::Dot(n, d) < 0)
    {
      std::cout << "Found a tetrahedra with bad winding." << std::endl;
      throw BoxClipTriangulateFailed();
    }
  }
}

//-----------------------------------------------------------------------------

static vtkSmartPointer<vtkUnstructuredGrid> BuildInput(
  int type, vtkIdType numcells, const vtkIdType* cells)
{
  vtkIdType i;

  VTK_CREATE(vtkUnstructuredGrid, input);

  // Randomly shuffle the points to possibly test various tessellations.
  // Make a map from original point orderings to new point orderings.
  std::vector<vtkIdType> idMap;
  std::vector<vtkIdType> idsLeft;

  for (i = 0; i < NumPoints; i++)
  {
    idsLeft.push_back(i);
  }

  while (!idsLeft.empty())
  {
    vtkIdType next = std::lround(vtkMath::Random(-0.49, idsLeft.size() - 0.51));
    std::vector<vtkIdType>::iterator nextp = idsLeft.begin() + next;
    idMap.push_back(*nextp);
    idsLeft.erase(nextp, nextp + 1);
  }

  // Build shuffled points.
  VTK_CREATE(vtkPoints, points);
  points->SetNumberOfPoints(NumPoints);
  for (i = 0; i < NumPoints; i++)
  {
    points->SetPoint(idMap[i], PointData + 3 * i);
  }
  input->SetPoints(points);

  // Add the cells with indices properly mapped.
  VTK_CREATE(vtkIdList, ptIds);
  const vtkIdType* c = cells;
  for (i = 0; i < numcells; i++)
  {
    vtkIdType npts = *c;
    c++;
    ptIds->Initialize();
    for (vtkIdType j = 0; j < npts; j++)
    {
      ptIds->InsertNextId(idMap[*c]);
      c++;
    }
    input->InsertNextCell(type, ptIds);
  }

  return input;
}

//-----------------------------------------------------------------------------

static void Check2DPrimitive(int type, vtkIdType numcells, const vtkIdType* cells)
{
  vtkSmartPointer<vtkUnstructuredGrid> input = BuildInput(type, numcells, cells);

  VTK_CREATE(vtkBoxClipDataSet, clipper);
  clipper->SetInputData(input);
  // Clip nothing.
  clipper->SetBoxClip(0.0, 2.0, 0.0, 1.0, 0.0, 1.0);
  clipper->Update();

  vtkUnstructuredGrid* output = clipper->GetOutput();

  if (output->GetNumberOfCells() < 1)
  {
    std::cout << "Output has no cells!" << std::endl;
    throw BoxClipTriangulateFailed();
  }

  // Check to make sure all the normals point in the z direction.
  vtkCellArray* outCells = output->GetCells();
  outCells->InitTraversal();
  vtkIdType npts;
  const vtkIdType* pts;
  while (outCells->GetNextCell(npts, pts))
  {
    if (npts != 3)
    {
      std::cout << "Got a primitive that is not a triangle!" << std::endl;
      throw BoxClipTriangulateFailed();
    }

    double n[3];
    vtkTriangle::ComputeNormal(output->GetPoints(), npts, pts, n);
    if ((n[0] > 0.1) || (n[1] > 0.1) || (n[2] < 0.9))
    {
      std::cout << "Primitive is facing the wrong way!" << std::endl;
      throw BoxClipTriangulateFailed();
    }
  }
}

//-----------------------------------------------------------------------------

static void Check3DPrimitive(
  int type, vtkIdType numcells, const vtkIdType* cells, vtkIdType numSurfacePolys)
{
  vtkSmartPointer<vtkUnstructuredGrid> input = BuildInput(type, numcells, cells);

  VTK_CREATE(vtkBoxClipDataSet, clipper);
  clipper->SetInputData(input);
  // Clip nothing.
  clipper->SetBoxClip(0.0, 2.0, 0.0, 1.0, 0.0, 1.0);
  clipper->Update();

  vtkUnstructuredGrid* output = clipper->GetOutput();

  if (output->GetNumberOfCells() < 1)
  {
    std::cout << "Output has no cells!" << std::endl;
    throw BoxClipTriangulateFailed();
  }

#if 0
  VTK_CREATE(vtkExtractEdges, edges);
  edges->SetInput(output);
  VTK_CREATE(vtkPolyDataMapper, mapper);
  mapper->SetInputConnection(0, edges->GetOutputPort(0));
  VTK_CREATE(vtkActor, actor);
  actor->SetMapper(mapper);
  VTK_CREATE(vtkRenderer, renderer);
  renderer->AddActor(actor);
  VTK_CREATE(vtkRenderWindow, renwin);
  renwin->AddRenderer(renderer);
  VTK_CREATE(vtkRenderWindowInteractor, iren);
  iren->SetRenderWindow(renwin);
  renwin->Render();
  iren->Start();
#endif

  CheckWinding(clipper);

  VTK_CREATE(vtkDataSetSurfaceFilter, surface);
  surface->SetInputConnection(clipper->GetOutputPort());
  surface->Update();

  if (surface->GetOutput()->GetNumberOfCells() != numSurfacePolys)
  {
    std::cout << "Expected " << numSurfacePolys << " triangles on the surface, got "
              << surface->GetOutput()->GetNumberOfCells() << std::endl;
    throw BoxClipTriangulateFailed();
  }
}

//-----------------------------------------------------------------------------

int BoxClipTriangulate(int, char*[])
{
  long seed = time(nullptr);
  std::cout << "Random seed = " << seed << std::endl;
  vtkMath::RandomSeed(seed);
  vtkMath::Random();
  vtkMath::Random();
  vtkMath::Random();

  try
  {
    std::cout << "Checking triangle strip." << std::endl;
    Check2DPrimitive(VTK_TRIANGLE_STRIP, NumTriStripCells, TriStripCells);

    std::cout << "Checking quadrilaterals." << std::endl;
    Check2DPrimitive(VTK_QUAD, NumQuadCells, QuadCells);

    std::cout << "Checking pixels." << std::endl;
    Check2DPrimitive(VTK_PIXEL, NumPixelCells, PixelCells);

    std::cout << "Checking polygons." << std::endl;
    Check2DPrimitive(VTK_POLYGON, NumPolyCells, PolyCells);

    std::cout << "Checking hexahedrons." << std::endl;
    Check3DPrimitive(VTK_HEXAHEDRON, NumHexCells, HexCells, NumExpectedHexSurfacePolys);

    std::cout << "Checking voxels." << std::endl;
    Check3DPrimitive(VTK_VOXEL, NumVoxelCells, VoxelCells, NumExpectedVoxelSurfacePolys);

    std::cout << "Checking wedges." << std::endl;
    Check3DPrimitive(VTK_WEDGE, NumWedgeCells, WedgeCells, NumExpectedWedgeSurfacePolys);

    std::cout << "Checking pyramids." << std::endl;
    Check3DPrimitive(VTK_PYRAMID, NumPyramidCells, PyramidCells, NumExpectedPyramidSurfacePolys);
  }
  catch (BoxClipTriangulateFailed)
  {
    return EXIT_FAILURE;
  }

  return EXIT_SUCCESS;
}
