//=========================================================================
//  Copyright (c) Kitware, Inc.
//  All rights reserved.
//  See LICENSE.txt 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.
//=========================================================================

#include "smtk/session/aeva/Registrar.h"
#include "smtk/session/aeva/Resource.h"
#include "smtk/session/aeva/Session.h"

#include "smtk/session/aeva/operators/Import.h"
#include "smtk/session/aeva/operators/ReconstructSurface.h"
#include "smtk/session/aeva/operators/SmoothSurface.h"

#include "smtk/attribute/Attribute.h"
#include "smtk/attribute/ComponentItem.h"
#include "smtk/attribute/DoubleItem.h"
#include "smtk/attribute/FileItem.h"
#include "smtk/attribute/IntItem.h"
#include "smtk/attribute/StringItem.h"

#include "smtk/model/EntityRef.h"
#include "smtk/model/Face.h"
#include "smtk/model/Model.h"

#include "smtk/resource/Manager.h"

#include "smtk/operation/Manager.h"

#include <vtkMath.h>
#include <vtkPoints.h>
#include <vtkPolyData.h>
#include <vtkTriangle.h>

namespace
{
std::string dataRoot = AEVA_DATA_DIR;

std::array<double, 2> angleStats(vtkPolyData* polydata)
{
  polydata->GetPolys()->InitTraversal();
  vtkSmartPointer<vtkIdList> idList = vtkSmartPointer<vtkIdList>::New();
  double p[3][3];
  double sum = 0.;
  double sum2 = 0.;
  int n = 0;
  for (vtkIdType i = 0; i < polydata->GetNumberOfCells(); i++)
  {
    vtkCell* cell = polydata->GetCell(i);
    if (vtkTriangle* triangle = vtkTriangle::SafeDownCast(cell))
    {
      triangle->GetPoints()->GetPoint(0, p[0]);
      triangle->GetPoints()->GetPoint(1, p[1]);
      triangle->GetPoints()->GetPoint(2, p[2]);
    }

    for (int j = 0; j < 3; j++)
    {
      std::array<double, 3> v1 = {
        p[(j + 1) % 3][0] - p[j][0], p[(j + 1) % 3][1] - p[j][1], p[(j + 1) % 3][2] - p[j][2]
      };
      std::array<double, 3> v2 = {
        p[(j + 2) % 3][0] - p[j][0], p[(j + 2) % 3][1] - p[j][1], p[(j + 2) % 3][2] - p[j][2]
      };

      double angle = vtkMath::AngleBetweenVectors(v1.data(), v2.data()) * 180. / vtkMath::Pi();
      sum += angle;
      sum2 += (angle * angle);
      ++n;
    }
  }

  double mean = sum / n;
  double variance = (sum2 - ((sum * sum) / n)) / n;
  return { mean, std::sqrt(variance) };
}
}

int TestSmoothSurfaceOp(int argc, char* argv[])
{
  (void)argc;
  (void)argv;

  // Create a resource manager
  smtk::resource::Manager::Ptr resourceManager = smtk::resource::Manager::create();

  // Register the aeva session to the resource manager
  {
    smtk::session::aeva::Registrar::registerTo(resourceManager);
  }

  // Create an operation manager
  smtk::operation::Manager::Ptr operationManager = smtk::operation::Manager::create();

  // Register the aeva session to the operation manager
  {
    smtk::session::aeva::Registrar::registerTo(operationManager);
  }

  // Register the resource manager to the operation manager (newly created
  // resources will be automatically registered to the resource manager).
  operationManager->registerResourceManager(resourceManager);

  smtk::model::Entity::Ptr model;

  {
    // Create an import operation
    smtk::session::aeva::Import::Ptr importOp =
      operationManager->create<smtk::session::aeva::Import>();
    if (!importOp)
    {
      std::cerr << "No import operation\n";
      return 1;
    }

    // Set the file path
    std::string importFilePath(dataRoot);
    importFilePath += "/vtk/oks003_TBB_AGS.vtk";
    importOp->parameters()->findFile("filename")->setValue(importFilePath);

    // Test the ability to operate
    if (!importOp->ableToOperate())
    {
      std::cerr << "Import operation unable to operate\n";
      return 1;
    }

    // Execute the operation
    smtk::operation::Operation::Result importOpResult = importOp->operate();

    // Retrieve the resulting model
    smtk::attribute::ComponentItemPtr componentItem =
      std::dynamic_pointer_cast<smtk::attribute::ComponentItem>(
        importOpResult->findComponent("created"));

    // Access the generated model
    model = std::dynamic_pointer_cast<smtk::model::Entity>(componentItem->value());

    // Test for success
    if (importOpResult->findInt("outcome")->value() !=
      static_cast<int>(smtk::operation::Operation::Outcome::SUCCEEDED))
    {
      std::cerr << "Import operation failed\n";
      return 1;
    }

    // Test model validity
    if (!model->referenceAs<smtk::model::Model>().isValid())
    {
      std::cerr << "Imported model is invalid\n";
      return 1;
    }
  }

  smtk::session::aeva::Session::Ptr session =
    std::static_pointer_cast<smtk::session::aeva::Resource>(model->modelResource())->session();

  // Access all of the model's faces
  smtk::model::EntityRefs faces =
    model->modelResource()->entitiesMatchingFlagsAs<smtk::model::EntityRefs>(smtk::model::FACE);
  if (faces.empty())
  {
    std::cerr << "No faces\n";
    return 1;
  }

  smtk::model::Face face = *faces.begin();

  vtkSmartPointer<vtkPolyData> facePD =
    vtkPolyData::SafeDownCast(session->findStorage(face.entity()));

  if (!facePD)
  {
    std::cerr << "No geometry for face\n";
    return 1;
  }

  if (facePD->GetNumberOfPoints() != 96340 || facePD->GetNumberOfCells() != 192676)
  {
    std::cerr << "Face geometry contains an unexpected number of points and/or cells\n";
    return 1;
  }

  {
    // Create an ReconstructSurface operation
    smtk::session::aeva::ReconstructSurface::Ptr reconstructSurface =
      operationManager->create<smtk::session::aeva::ReconstructSurface>();
    if (!reconstructSurface)
    {
      std::cerr << "No reconstruct surface operation\n";
      return 1;
    }

    // Set the input face
    reconstructSurface->parameters()->associate(face.component());

    // Set the voxel grid dimensions
    reconstructSurface->parameters()->findInt("dimensions")->setValue(0, 100);
    reconstructSurface->parameters()->findInt("dimensions")->setValue(1, 100);
    reconstructSurface->parameters()->findInt("dimensions")->setValue(2, 100);

    // Set the radius of influence
    reconstructSurface->parameters()->findDouble("radius")->setValue(5.);

    // Test the ability to operate
    if (!reconstructSurface->ableToOperate())
    {
      std::cerr << "Reconstruct surface operation unable to operate\n";
      return 1;
    }

    // Execute the operation
    smtk::operation::Operation::Result reconstructSurfaceResult = reconstructSurface->operate();

    // Retrieve the resulting face
    smtk::attribute::ComponentItemPtr componentItem =
      std::dynamic_pointer_cast<smtk::attribute::ComponentItem>(
        reconstructSurfaceResult->findComponent("created"));

    // Access the generated model
    face = std::dynamic_pointer_cast<smtk::model::Entity>(componentItem->value())
             ->referenceAs<smtk::model::Face>();

    // Test for success
    if (reconstructSurfaceResult->findInt("outcome")->value() !=
      static_cast<int>(smtk::operation::Operation::Outcome::SUCCEEDED))
    {
      std::cerr << "Reconstruct surface operation failed\n";
      return 1;
    }

    // Test model validity
    if (!face.isValid())
    {
      std::cerr << "Reconstructed face is invalid\n";
      return 1;
    }
  }

  facePD = vtkPolyData::SafeDownCast(session->findStorage(face.entity()));

  if (!facePD)
  {
    std::cerr << "No geometry for face\n";
    return 1;
  }

  if (facePD->GetNumberOfPoints() != 62748 || facePD->GetNumberOfCells() != 76852)
  {
    std::cerr << "Face geometry contains an unexpected number of points and/or cells\n";
    return 1;
  }

  double initialSD;
  {
    std::array<double, 2> stats = angleStats(facePD);
    initialSD = stats[1];
  }

  {
    // Create an SmoothSurface operation
    smtk::session::aeva::SmoothSurface::Ptr smoothSurface =
      operationManager->create<smtk::session::aeva::SmoothSurface>();
    if (!smoothSurface)
    {
      std::cerr << "No smooth surface operation\n";
      return 1;
    }

    // Set the input face
    smoothSurface->parameters()->associate(face.component());

    // Test the ability to operate
    if (!smoothSurface->ableToOperate())
    {
      std::cerr << "Smooth surface operation unable to operate\n";
      return 1;
    }

    // Execute the operation
    smtk::operation::Operation::Result smoothSurfaceResult = smoothSurface->operate();

    // Retrieve the resulting face
    smtk::attribute::ComponentItemPtr componentItem =
      std::dynamic_pointer_cast<smtk::attribute::ComponentItem>(
        smoothSurfaceResult->findComponent("modified"));

    // Access the generated model
    face = std::dynamic_pointer_cast<smtk::model::Entity>(componentItem->value())
             ->referenceAs<smtk::model::Face>();

    // Test for success
    if (smoothSurfaceResult->findInt("outcome")->value() !=
      static_cast<int>(smtk::operation::Operation::Outcome::SUCCEEDED))
    {
      std::cerr << "Smooth surface operation failed\n";
      return 1;
    }

    // Test model validity
    if (!face.isValid())
    {
      std::cerr << "Smoothed face is invalid\n";
      return 1;
    }
  }

  facePD = vtkPolyData::SafeDownCast(session->findStorage(face.entity()));

  if (!facePD)
  {
    std::cerr << "No geometry for face\n";
    return 1;
  }

  if (facePD->GetNumberOfPoints() != 62748 || facePD->GetNumberOfCells() != 76852)
  {
    std::cerr << "Face geometry contains an unexpected number of points and/or cells\n";
    return 1;
  }

  double finalSD;
  {
    std::array<double, 2> stats = angleStats(facePD);
    finalSD = stats[1];
  }

  if (finalSD >= initialSD)
  {
    std::cerr
      << "Smoothing operation did not lower the standard deviation of the face triangles' angles\n";
    return 1;
  }

  return 0;
}
