#include "vtkLimitFibersToROI.h"

#include <vtkPointData.h>
#include <vtkDataArray.h>
#include <vtkPolyData.h>
#include <vtkObjectFactory.h>
#include <vtkInformationVector.h>
#include <vtkInformation.h>
#include <vtkCellArray.h>

#ifdef WIN32
extern int rint(double a);
#endif

vtkCxxRevisionMacro(vtkLimitFibersToROI, "$Revision: 1.12 $");
vtkStandardNewMacro(vtkLimitFibersToROI);

vtkLimitFibersToROI::vtkLimitFibersToROI()
{
  MaskImage = 0;
}

vtkLimitFibersToROI::~vtkLimitFibersToROI()
{
  if( MaskImage )
    MaskImage->Delete();
}


int vtkLimitFibersToROI::RequestData (vtkInformation *vtkNotUsed(request),
                                       vtkInformationVector **inputVector,
                                       vtkInformationVector *outputVector)
{

  // get the info objects
  vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
  vtkInformation *outInfo = outputVector->GetInformationObject(0);

  // get the input and ouptut
  vtkPolyData *input = vtkPolyData::SafeDownCast(inInfo->Get(vtkDataObject::DATA_OBJECT()));
  vtkPolyData *output = vtkPolyData::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT()));
  output->Initialize();
  output->Allocate();
  
  vtkPoints* points = input->GetPoints();
  output->SetPoints ( points );
  vtkDataArray* scalars = input->GetPointData()->GetScalars();
  output->GetPointData()->SetScalars ( scalars );
  
  vtkCellArray* lines = input->GetLines();
  if( lines==0 )
    return 0;
  
  if (MaskImage == 0)
  {
    output->SetLines (lines);
    return 1;
  }

  //store informations
  int* dim = MaskImage->GetDimensions();
  double* origin = MaskImage->GetOrigin();
  double* spacing = MaskImage->GetSpacing();
  std::vector<int> Labels;
  unsigned char* MaskValues= (unsigned char*)MaskImage->GetScalarPointer();
  
  for (int i=0;i<MaskImage->GetNumberOfPoints();i++)
  {
    int val = (int)(*MaskValues);    
    if (val>0)
    {
      bool IsInserted = false;
      for (unsigned int j=0;j<Labels.size();j++)
        if (Labels[j]==val)
          IsInserted=true;
      if (!IsInserted)
        Labels.push_back (val);
    }
    ++MaskValues;
  }
  std::cout<<"Number Of Valid ROIs: "<<Labels.size()<<std::endl;

  //reset the Scalar Pointer
  MaskValues= (unsigned char*)MaskImage->GetScalarPointer();
  
  //Algorithm
  
  // 2 Temporary polydata to carry the information
  vtkPolyData* Aux1 = vtkPolyData::New();
  vtkPolyData* Aux2 = vtkPolyData::New();

  Aux1->Initialize();
  Aux1->Allocate();
  Aux2->Initialize();
  Aux2->Allocate();
  Aux1->SetLines (lines);
  char temp[512];
  
  // for all the ROIs
  for (unsigned int i=0;i<Labels.size();i++)
  {
    vtkCellArray* auxLines = Aux1->GetLines();
    auxLines->InitTraversal();
        
    int  npts=0;
    int* ptids=0;
    int test = auxLines->GetNextCell (npts, ptids);

    sprintf (temp,"%d / %d", i+1, (int)Labels.size());
    //this->SetProgressText (temp);
    this->UpdateProgress (0.0);
    int iterProgress = 0;
    int maxProgress = Aux1->GetLines()->GetNumberOfCells();
    int step=maxProgress/100;
    
    // for All the Fibers 
    while ( test )
    {
      bool found = false;
      
      // for All the points of the fiber
      for (int k=0; k<npts; k++)
      {
        double pt[3];
        points->GetPoint (ptids[k], pt);
        
        int c[3];
        c[0] = (int)( rint( (pt[0]-origin[0])/spacing[0] ) );
        c[1] = (int)( rint( (pt[1]-origin[1])/spacing[1] ) );
        c[2] = (int)( rint( (pt[2]-origin[2])/spacing[2] ) );
        
        if( c[0]>=0 && c[0]<dim[0] &&
            c[1]>=0 && c[1]<dim[1] &&
            c[2]>=0 && c[2]<dim[2] )
        {
          int index = c[0] + c[1]*dim[0] + c[2]*dim[0]*dim[1];
          MaskValues+=index;
          int maskVal = (int)(*MaskValues);
          MaskValues-=index;
          
          if( maskVal==Labels[i] )
          {
            found=true;
            break;
          }
        }
      }
      
      if (found)
        Aux2->InsertNextCell(VTK_POLY_LINE, npts, ptids);
      
      test = auxLines->GetNextCell (npts, ptids);

      iterProgress++;
      if( step>0 )
        if ( !(iterProgress%step) )
          this->UpdateProgress (double(iterProgress)/double(maxProgress));
    }

    auxLines = Aux2->GetLines();
    Aux1->SetLines ( auxLines );

    Aux2->Initialize();
    Aux2->Allocate();     
  }
  
  output->SetLines ( Aux1->GetLines() );
  
  Aux1->Delete();
  Aux2->Delete();
    
  this->UpdateProgress (1.0);
  

  return 1;
  
}