#ifndef MESHSTLWRITER_H_HEADER_INCLUDED
#define MESHSTLWRITER_H_HEADER_INCLUDED
#include "itkMeshFileWriter.h"
#include <vnl/vnl_vector_fixed.h>
#include <vnl/vnl_cross.h>

namespace itk {

/** \class MeshSTLWriter
  *   \brief Stores a mesh in binary STL (stereolitography) format.
  *
  * Only triangular faces are stored. 
  * Use SetFilePrefix() to specify the filename, which will be extended with 
  * the stl-extensions.
  *
  * \author Tobias Heimann. Division Medical and Biological Informatics, 
  *         German Cancer Research Center, Heidelberg, Germany.
  */
template <class TMeshType>
class MeshSTLWriter : public MeshFileWriter<TMeshType>
{
  public:

    typedef MeshSTLWriter<TMeshType> Self;
    typedef MeshFileWriter<TMeshType> Superclass;
    typedef SmartPointer < Self > Pointer;
    typedef SmartPointer < const Self > ConstPointer;

    itkTypeMacro( MeshSTLWriter<TMeshType>, MeshFileWriter<TMeshType> );
    itkNewMacro( Self );  


    virtual void GenerateData()
    {
      if (strlen( this->m_filePrefix ) == 0) {
        itkExceptionMacro(<<"GenerateData:Please specify a file prefix!");
        return;
      }
      char fileName[1024];
      typedef typename TMeshType::CellsContainer::ConstIterator CellIterator;
      typedef typename TMeshType::CellType CellType;
      typename TMeshType::Pointer mesh = dynamic_cast<TMeshType*> (this->GetInput());
      
      sprintf( fileName, "%s.stl", this->m_filePrefix );
      FILE *file = fopen( fileName, "wb" );
      if (file) {
        // count triangles in mesh
        unsigned int numTriangles = 0;
        CellIterator cellIterator = mesh->GetCells()->Begin();
        CellIterator cellEnd      = mesh->GetCells()->End();
        while( cellIterator != cellEnd ) {
          CellType *cell = cellIterator.Value();
          if (cell->GetType() == CellType::TRIANGLE_CELL) numTriangles++;
          cellIterator++;
        }
        // write header
        char header[81];
        strcpy( header, "STL generated by itk::MeshSTLWriter class. Author: t.heimann@dkfz.de." );
        fwrite( (void*)&(header[0]), 1, 80, file );
        fwrite( &numTriangles, 4, 1, file );
        // write poly data
        float vertex[3][3];
        vnl_vector_fixed<float, 3> v1, v2, normal;
        short stlSeperator = 0;
        typename TMeshType::PointType meshPoint;
        cellIterator = mesh->GetCells()->Begin();
        while( cellIterator != cellEnd ) {
          CellType *cell = cellIterator.Value();
          if (cell->GetType() == CellType::TRIANGLE_CELL) {
            typename CellType::PointIdConstIterator pntId = cell->GetPointIds();
            for (int i=0; i<3; i++) {
              mesh->GetPoint( pntId[i], &meshPoint );
              for (int d=0; d<3; d++) vertex[i][d] = meshPoint[d];
            }
            for (int d=0; d<3; d++) {
              v1[d] = vertex[1][d] - vertex[0][d];  
              v2[d] = vertex[2][d] - vertex[0][d];  
            }
            normal = vnl_cross_3d( v1, v2 );
            normal /= normal.magnitude();
    
            fwrite( (void*)normal.data_block(), sizeof(float), 3, file );
            fwrite( (void*)&(vertex[0][0]), sizeof(float), 9, file );
            fwrite( (void*)&(stlSeperator), 2, 1, file );
          }
          cellIterator++;
        }
      }
      fclose( file );
    }

   
protected:
  MeshSTLWriter() {};
  virtual ~MeshSTLWriter() {};

};

} // namespace itk
#endif /* MESHSTLWRITER_H_HEADER_INCLUDED */