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

  Program:   Visualization Toolkit
  Module:    vtkGeneralizedSurfaceNets3D.h

  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.

=========================================================================*/
/**
 * @class   vtkGeneralizedSurfaceNets3D
 * @brief   create a surface net from a unorganized set of segmented (i.e., labeled) points
 *
 * vtkGeneralizedSurfaceNets3D is a filter that constructs a surface net from
 * a labeled / segmentedlist of input points. The points are presumed to lie
 * within 3D-space. These points may be represented by any dataset of type
 * vtkPointSet and subclasses. The output of the filter is a complex of
 * convex polygons represented by a vtkPolyData. Additionally the output
 * contains cell data consisting of a 2-component tuples which record the
 * regions on either side of the polygonal faces composing the surface
 * net. The algorithm uses a novel 3D Voronoi tessellation algorithm, and
 * extracts surface net faces between Voronoi hulls which lie in specified
 * separate regions.
 *
 * Besides the list of input points, the filter requires an input region ids
 * array which labels the points as belonging to different regions. This
 * point data region ids array should be of integer type, with region id
 * values>=0. Note that a region id < 0 indicates that the associated point
 * is "outside"; consequently the point will not produce output but will
 * affect neighboring points in terms of producing boundary faces.
 *
 * The surface net algorithm can also (optionally) smooth the output
 * polygonal surface. To be faithful to the original algorithm, a
 * vtkConstrainedSmoothingFilter is used; however other smoothing algorithms
 * such as vtkWindowedSincPolyDataFilter may be used (by disabling smoothing,
 * enabling point merging, and passing the output of the filter to subseguent
 * smoothing filter).
 *
 * There are two common use cases when using this filter. The first case
 * simply produces output surface net faces for the purposes of
 * visualization. In this case the resulting surfaces are not watertight and
 * cannot be smoothed (so-called meshless complex of polygons). (Note that
 * non-smoothed surface nets tend to be choppy depending on the input point
 * cloud resolution.) The second case produces connected, watertight surface
 * meshes which can be smoothed. Note that this second case requires a fair
 * amount of work to merge nearly coincident points to produce the watertight
 * surfaces. (Note: a built-in topologically-based point merging process is
 * used. Users can disable the built in point merging process, and use
 * subsequent filters like vtkStaticCleanPolyData to merge coincident points,
 * remove degenerate face primitives, etc, and otherwise process the surfaces
 * with smoothing etc. vtkStaticCleanPolyData uses a proximal geometric point
 * merging process requiring a tolerance, this can cause problems in some
 * cases.)
 *
 * Finally, another important option to the filter is that surfaces
 * corresponding to the domain boundary can be generated. In some cases it is
 * useful to provide the boundary as context to the surface net
 * contours. (The domain boundary is determined from Voronoi edges that
 * connect to the domain edges, or connect to points with region ids < 0).
 *
 * See the following reference for more details about surface nets:
 * W. Schroeder, S. Tsalikis, M. Halle, S. Frisken. A High-Performance
 * SurfaceNets Discrete Isocontouring Algorithm. arXiv:2401.14906. 2024.
 * (http://arxiv.org/abs/2401.14906).
 *
 * The Surface Nets algorithm was first proposed by Sarah Frisken.  Two
 * important papers include the description of surface nets for binary
 * objects (i.e., extracting just one segmented object from a volume) and
 * multi-label (multiple object extraction).
 *
 * S. Frisken (Gibson), “Constrained Elastic SurfaceNets: Generating Smooth
 * Surfaces from Binary Segmented Data”, Proc. MICCAI, 1998, pp. 888-898.
 *
 * S. Frisken, “SurfaceNets for Multi-Label Segmentations with Preservation
 * of Sharp Boundaries”, J. Computer Graphics Techniques, 2022.
 *
 * These techniques mentioned above are specialized to input 3D labeled (or
 * segmented) volumes. This filter implementes a generalized version of
 * surface nets for labeled, unorganized point clouds.
 *
 * @warning
 * Coincident input points are discarded. The Voronoi tessellation requires
 * unique input points.
 *
 * @warning
 * This class has been threaded with vtkSMPTools. Using TBB or other
 * non-sequential type (set in the CMake variable
 * VTK_SMP_IMPLEMENTATION_TYPE) may improve performance significantly.
 *
 * @sa
 * vtkSurfaceNets2D vtkSurfaceNets3D vtkVoronoi3D
 * vtkConstrainedSmoothingFilter vtkWindowedSincPolyDataFilter
 */

#ifndef vtkGeneralizedSurfaceNets3D_h
#define vtkGeneralizedSurfaceNets3D_h

#include "vtkConstrainedSmoothingFilter.h" // Perform surface smoothing
#include "vtkContourValues.h" // Manage countour values
#include "vtkFiltersCoreModule.h" // For export macro
#include "vtkStaticPointLocator.h" // For point locator
#include "vtkPolyDataAlgorithm.h"


class VTKFILTERSCORE_EXPORT vtkGeneralizedSurfaceNets3D : public vtkPolyDataAlgorithm
{
public:
  ///@{
  /**
   * Standard methods for instantiation, type information, and printing.
   */
  static vtkGeneralizedSurfaceNets3D* New();
  vtkTypeMacro(vtkGeneralizedSurfaceNets3D, vtkPolyDataAlgorithm);
  void PrintSelf(ostream& os, vtkIndent indent) override;
  ///@}

  //---------------The following are methods used to set/get and generate
  //---------------contour labels.
   ///@{
  /**
   * Set a particular label value at label number i. The index i ranges
   * between (0 <= i < NumberOfLabels). (Note: while labels values are
   * expressed as doubles, the underlying scalar data may be a different
   * type. During execution the label values are cast to the type of the
   * scalar data.)  Note the use of "Value" and "Label" when specifying
   * regions to extract. The use of "Value" is consistent with other VTK
   * continuous-scalar field isocontouring algorithms; however the term
   * "Label" is more consistent with label maps.  Warning: make sure that the
   * label value >= 0 as any label < 0 is considered a background, i.e.,
   * outside, label.
   */
  void SetValue(int i, double value) { this->Labels->SetValue(i, value); }
  void SetLabel(int i, double value) { this->Labels->SetValue(i, value); }
  ///@}

  ///@{
  /**
   * Get the ith label value.
   */
  double GetValue(int i) { return this->Labels->GetValue(i); }
  double GetLabel(int i) { return this->Labels->GetValue(i); }
  ///@}

  ///@{
  /**
   * Get a pointer to an array of labels. There will be
   * GetNumberOfLabels() values in the list.
   */
  double* GetValues() { return this->Labels->GetValues(); }
  double* GetLabels() { return this->Labels->GetValues(); }
  ///@}

  ///@{
  /**
   * Fill a supplied list with label values. There will be
   * GetNumberOfLabels() values in the list. Make sure you allocate enough
   * memory to hold the list.
   */
  void GetValues(double* contourValues) { this->Labels->GetValues(contourValues); }
  void GetLabels(double* contourValues) { this->Labels->GetValues(contourValues); }
  ///@}

  ///@{
  /**
   * Set the number of labels to place into the list. You only really need to
   * use this method to reduce list size. The method SetValue() will
   * automatically increase list size as needed. Note that for consistency
   * with other isocountoring-related algorithms, some methods use
   * "Labels" and "Contours" interchangeably.
   */
  void SetNumberOfLabels(int number) { this->Labels->SetNumberOfContours(number); }
  void SetNumberOfContours(int number) { this->Labels->SetNumberOfContours(number); }
  ///@}

  ///@{
  /**
   * Get the number of labels in the list of label values.
   */
  vtkIdType GetNumberOfLabels() { return this->Labels->GetNumberOfContours(); }
  vtkIdType GetNumberOfContours() { return this->Labels->GetNumberOfContours(); }
  ///@}

  ///@{
  /**
   * Generate numLabels equally spaced labels between the specified
   * range. The labels will include the min/max range values.
   */
  void GenerateLabels(int numLabels, double range[2])
  {
    this->Labels->GenerateValues(numLabels, range);
  }
  void GenerateValues(int numContours, double range[2])
  {
    this->Labels->GenerateValues(numContours, range);
  }
  void GenerateLabels(int numLabels, double rangeStart, double rangeEnd)
  {
    this->Labels->GenerateValues(numLabels, rangeStart, rangeEnd);
  }
  void GenerateValues(int numContours, double rangeStart, double rangeEnd)
  {
    this->Labels->GenerateValues(numContours, rangeStart, rangeEnd);
  }
  ///@}

  ///@{
  /**
   * This value specifies the label value to use when indicating that a region
   * is outside. That is, the output 2-tuple cell data array indicates which
   * segmented region is on either side of it. To indicate that one side is
   * on the boundary, the BackgroundLabel value is used. The BackgroundLabel
   * must be <0. By default the background label is -100.
   */
  vtkSetClampMacro(BackgroundLabel, int, VTK_INT_MIN, -1);
  vtkGetMacro(BackgroundLabel, int);
  ///@}

  //---------------Done defining label-related methods.
  ///@{
  /**
   * Specify whether to cap the surface net along the boundary.  By default
   * this is off.
   */
  vtkGetMacro(BoundaryCapping, vtkTypeBool);
  vtkSetMacro(BoundaryCapping, vtkTypeBool);
  vtkBooleanMacro(BoundaryCapping, vtkTypeBool);
  ///@}

  ///@{
  /**
   * Specify whether to merge nearly concident points in order to produce
   * watertight output surfaces. Enabling merging is necessary to perform
   * smoothing. However it does require significant time to compute. By
   * default this is on.
   */
  vtkGetMacro(MergePoints, vtkTypeBool);
  vtkSetMacro(MergePoints, vtkTypeBool);
  vtkBooleanMacro(MergePoints, vtkTypeBool);
  ///@}

  ///@{
  /**
   * Indicate whether smoothing should be enabled. By default, after the
   * surface net is extracted, smoothing occurs using the built-in smoother,
   * and MergePoints is enabled.  To disable smoothing, invoke
   * SmoothingOff(). (Note: disabling smoothing can be useful if simply
   * visualizing the initial surface net, or if you desire to use a different
   * smoother.)
   */
  vtkSetMacro(Smoothing, vtkTypeBool);
  vtkGetMacro(Smoothing, vtkTypeBool);
  vtkBooleanMacro(Smoothing, vtkTypeBool);
  ///@}

  ///@{
  /**
   * Convenience methods that delegate to the internal smoothing filter
   * follow below. See the documentation for vtkConstrainedSmoothingAlgorithm
   * for more information.
   */
  void SetNumberOfIterations(int n) { this->Smoother->SetNumberOfIterations(n); }
  int GetNumberOfIterations() { return this->Smoother->GetNumberOfIterations(); }
  void SetRelaxationFactor(double f) { this->Smoother->SetRelaxationFactor(f); }
  double GetRelaxationFactor() { return this->Smoother->GetRelaxationFactor(); }
  void SetConstraintDistance(double d) { this->Smoother->SetConstraintDistance(d); }
  double GetConstraintDistance() { return this->Smoother->GetConstraintDistance(); }
  ///@}

  ///@{
  /**
   * Get the instance of vtkConstrainedSmoothingFilter used to smooth the
   * extracted surface net. To control smoothing, access this instance and
   * specify its parameters such as number of smoothing iterations and
   * constraint distance. If you wish to disable smoothing, set
   * SmoothingOff().
   */
  vtkGetSmartPointerMacro(Smoother, vtkConstrainedSmoothingFilter);
  ///@}

  ///@{
  /**
   * Specify a padding for the bounding box of the input points. A >0 padding
   * is necessary in order to create valid Voronoi hulls on the boundary of
   * the tessellation. The padding is specified as a fraction of the diagonal
   * length of the bounding box of the points.
   */
  vtkSetClampMacro(Padding, double, 0.001, 0.25);
  vtkGetMacro(Padding, double);
  ///@}

  ///@{
  /**
   * Retrieve the internal locator to manually configure it, for example
   * specifying the number of points per bucket. This method is generally
   * used for debugging or testing purposes.
   */
  vtkStaticPointLocator* GetLocator() { return this->Locator; }
  ///@}

  ///@{
  /**
   * These methods are for debugging or instructional purposes. When the
   * point of interest is specified (i.e., set to a non-negative number) then
   * the algorithm will process this single point (whose id is the
   * PointOfInterest). When PointsOfInterest is specified through a supplied
   * vtkIdTypeArray (this is in addition to the PointOfInterest), then only
   * those hulls in the PointOfInterest + PointsOfInterestArray will be
   * produced. The maximum number of clips (the MaximumNumberOfHullClips) can
   * be specified. If MaximumNumberOfHullClips=0, then the initial tile
   * (single point within the bounding box) is produced; if =1 then the split
   * with the closest point is produced; and so on. By default the
   * PointOfInterest is set to (-1), and the number of clips is unlimited
   * (i.e., MaximumNumberOfHullClips=VTK_ID_MAX and therefore automatically
   * limited by the algorithm).
   */
  vtkSetClampMacro(PointOfInterest, vtkIdType, -1, VTK_ID_MAX);
  vtkGetMacro(PointOfInterest, vtkIdType);
  vtkSetObjectMacro(PointsOfInterest,vtkIdTypeArray);
  vtkGetObjectMacro(PointsOfInterest,vtkIdTypeArray);
  vtkSetClampMacro(MaximumNumberOfHullClips, vtkIdType, 1, VTK_ID_MAX);
  vtkGetMacro(MaximumNumberOfHullClips, vtkIdType);
  ///@}

  ///@{
  /**
   * Specify the number of input generating points in a batch, where a batch
   * defines a contiguous subset of the input points operated on during
   * threaded execution. Generally this is only used for debugging or
   * performance studies (since batch size affects the thread workload).
   *
   * Default is 1000.
   */
  vtkSetClampMacro(BatchSize, unsigned int, 1, VTK_INT_MAX);
  vtkGetMacro(BatchSize, unsigned int);
  ///@}

  /**
   *  Return the number of threads actually used during execution. This is
   *  valid only after algorithm execution.
   */
  int GetNumberOfThreadsUsed() { return this->NumberOfThreadsUsed; }

 /**
   * The modified time is also a function of the built in locator, smoothing
   * filter, and label values.
   */
   vtkMTimeType GetMTime() override;

protected:
  vtkGeneralizedSurfaceNets3D();
  ~vtkGeneralizedSurfaceNets3D() override = default;

  // Support the contouring operation by defining labels.
  vtkSmartPointer<vtkContourValues> Labels;
  int BackgroundLabel;

  // Algorithm control
  vtkTypeBool BoundaryCapping; // produce boundary surfaces or not
  vtkTypeBool MergePoints; // merge near coincident points or not
  vtkTypeBool Smoothing; // enable built-in smoothing process

  // Related to smoothing control
  vtkSmartPointer<vtkConstrainedSmoothingFilter> Smoother; // built in smoother
  vtkSmartPointer<vtkStaticPointLocator> Locator; // locator for finding proximal points

  // Related to internal Voronoi methods
  double Padding; // amount to pad out input points bounding box
  vtkIdType PointOfInterest; // specify a single input point to process
  vtkSmartPointer<vtkIdTypeArray> PointsOfInterest; // list of points of interest
  vtkIdType MaximumNumberOfHullClips; // limit the number of hull clips
  unsigned int BatchSize; // process data in batches of specified size
  int NumberOfThreadsUsed; // report on the number of threads used during processing

  // Satisfy pipeline-related API
  int RequestData(vtkInformation*, vtkInformationVector**, vtkInformationVector*) override;
  int FillInputPortInformation(int port, vtkInformation* info) override;

private:
  vtkGeneralizedSurfaceNets3D(const vtkGeneralizedSurfaceNets3D&) = delete;
  void operator=(const vtkGeneralizedSurfaceNets3D&) = delete;
};

#endif
