/*---------------------------------------------------------------------------*\

    geophase: A C++ template library for geometrical operations on non-convex
    polyhedra with non-planar faces
     
    Copyright (C) 2019 Tomislav Maric, TU Darmstadt 

    This program is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program.  If not, see <http://www.gnu.org/licenses/>.

Description
    Legacy VTK POLYDATA ASCII output. 

Author
    Tomislav Maric maric@mma.tu-darmstadt.de

\*---------------------------------------------------------------------------*/


#ifndef WriteVtkPolyData_HPP
#define WriteVtkPolyData_HPP

#include <fstream>
#include <sstream>
#include <iomanip>

#include "tag.hpp"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace geophase { 

template<typename Stream>
class VtkPolyDataOStream 
{
    std::size_t precision_; 
    Stream pointStream_; 
    //Stream vertexStream_; 
    Stream polygonStream_; 
    Stream fieldStream_; 
    std::ofstream outFile_; 

    unsigned long nPoints_ = 0; 
    unsigned long nPolygons_ = 0;
    unsigned long polygonDimension_ = 0; 

    public: 

        VtkPolyDataOStream() = default; 

        template<typename FileName> 
        VtkPolyDataOStream(
            FileName const& fileName, 
            const char* dataName = "polyData", 
            std::size_t precision = 17
        )
        :
            precision_(precision),
            pointStream_(),
            polygonStream_(),
            fieldStream_(),
            outFile_()
        {
            outFile_.open(fileName); 
            writeHeader(outFile_, dataName); 
        }

        ~VtkPolyDataOStream()
        {
            close(outFile_);  
        }

        template<typename OutputStream, typename Name>
        void writeHeader(OutputStream& os, Name const & name)
        {
            os << "# vtk DataFile Version 2.0\n" << name << "\nASCII\nDATASET POLYDATA\n"; 
        } 

        template<typename OutputStream>
        void close(OutputStream& os)
        {
            os  << "POINTS " << nPoints_ << " float\n" 
                << pointStream_.str() << "\n" 
                //<< "VERTICES " << nPoints_ << " " << nPoints_ * 2 << "\n" 
                //<< vertexStream_.str() << "\n"
                << "POLYGONS " << nPolygons_ << " " << polygonDimension_ 
                << "\n" << polygonStream_.str() << "\n";
            os.close(); 
        }

        void close()
        {
            close(outFile_); 
        }

        void setPrecision(std::size_t precision)
        {
            precision_ = precision;
        }

        template<typename Geometry>
        void streamGeometry(Geometry const& geo)
        {
            using tag = typename tag<Geometry>::type; 
            streamGeometry(geo, tag());
        }

        template<typename Point>
        void streamGeometry(Point const& point, vector_tag)
        {
            pointStream_ << std::setprecision(precision_) 
                << point[0] << " " << point[1] 
                << " " << point[2] << "\n"; 
            
            //vertexStream_ << 1 << " " << nPoints_ << "\n"; 

            ++nPoints_; 
        }

        //template<typename PointSequence>
        //void streamGeometry(PointSequence const& points, sequence_tag)
        //{
            //for (const auto& point : points)
            //{
                //pointStream_ << std::setprecision(precision_) 
                    //<< point[0] << " " << point[1] 
                    //<< " " << point[2] << "\n"; 

                ////vertexStream_ << 1 << " " << nPoints_ << "\n"; 
                //++nPoints_; 
            //}
        //}

        template<typename Polygon>
        void streamGeometry(Polygon const& polygon, sequence_polygon_tag)
        {
            std::size_t polygonSize = polygon.size(); 
            auto polygonPointIndex = nPoints_;  
            polygonStream_ << polygonSize << " "; 
            for (const auto& point : polygon)
            {
                streamGeometry(point); 
                polygonStream_ << polygonPointIndex << " ";  
                ++polygonPointIndex; 
            }
            polygonStream_ << "\n"; 
            ++nPolygons_; 
            polygonDimension_ += polygonSize + 1;  
        }

        template<typename PolygonCollection>
        void streamGeometry(
            PolygonCollection const& polyhedron, 
            sequence_polyhedron_tag
        )
        {
            for (const auto& polygon : polyhedron)
                streamGeometry(polygon); 
        }

        template<typename GeometrySequence>
        void streamGeometry(
            GeometrySequence const& seq, 
            sequence_tag 
        )
        {
            for (const auto& geo : seq)
                streamGeometry(geo); 
        }

        //template<typename Box>
        //void streamGeometry(Box const& box, aabbox_tag)
        //{
            //streamGeometry(build<pointVectorVector>(box)); 
        //}

        //template<typename Triangle>
        //void streamGeometry(Triangle const& triangle, triangle_tag)
        //{
            //streamGeometry(triangle, point_collection_tag()); 
        //}

        //template<typename Tetrahedron>
        //void streamGeometry(Tetrahedron const& tet, tetrahedron_tag)
        //{
            //for(char I = 0; I < 4; ++I)
            //{
                //streamGeometry(tet[I], point_collection_tag()); 
            //} 
        //}

        //template<typename Triangulation>
        //void streamGeometry(Triangulation& tri, triangle_sequence_tag)
        //{
            //for(const auto& triangle : tri)
            //{
                //streamGeometry(triangle); 
            //} 
        //}

        //template<typename Triangulation>
        //void streamGeometry(Triangulation& tri, tetrahedron_sequence_tag)
        //{
            //for(const auto& tetrahedron : tri)
            //{
                //streamGeometry(tetrahedron); 
            //} 
        //}

        //template<typename TriangulationIntersection>
        //void streamGeometry(TriangulationIntersection& intersection, triangulation_intersection_tag)
        //{
            //for(const auto& tetrahedron : intersection.tetrahedra())
            //{
                //streamGeometry(tetrahedron); 
            //} 

            //for(const auto& polyhedron : intersection.polyhedra())
            //{
                //streamGeometry(polyhedron); 
            //} 
        //}

        //template<typename MixedPolyhedron>
        //void streamGeometry(MixedPolyhedron& mixedPolyhedron, mixed_polyhedron_tag)
        //{
            //for(const auto& polygon : mixedPolyhedron.polygons())
            //{
                //streamGeometry(polygon); 
            //} 

            //for(const auto& triangle : mixedPolyhedron.triangles())
            //{
                //streamGeometry(triangle); 
            //} 
        //}

        //template<typename Geometry>
        //void streamGeometry(Geometry const& geo, pair_tag)
        //{
            //streamGeometry(geo.second); 
        //}

        //template<typename Geometry>
        //void streamGeometry(Geometry const& interface, geometrical_interface_tag)
        //{
            //for (const auto& element : interface)
            //{
                //// TODO: remove this? Enforce correct polyhedron sizes? 
                //const auto& polyhedron = element.polyhedron(); 
                //if (polyhedron.size() > 0)
                //{
                    //streamGeometry(element.polygon()); 
                //}
            //}
        //}


        //template<typename ScalarField, typename String> 
        //void streamField(const ScalarField& sf, String name) 
        //{
            //++nFields;  
            //fieldStream_ << name << " " << 1 << " " << sf.size() << " float \n";
            //for (const auto& x : sf)
            //{
                //fieldStream_ << x << " "; 
            //}
            //fieldStream_ << "\n"; 
        //} 
};

typedef VtkPolyDataOStream<std::stringstream> vtkPolyDataOStream;  

template<typename Geometry> 
vtkPolyDataOStream& operator<<(
    vtkPolyDataOStream& polyDataStream, 
    Geometry const & geo
)
{
    polyDataStream.streamGeometry(geo); 
    return polyDataStream; 
}

template <typename Geometry>
void write_vtk(Geometry const& geom, std::string const& fileName)
{
    vtkPolyDataOStream ps(fileName);  
    ps.streamGeometry(geom); 
}

template<typename P, typename ... Q>
void stream_labels(std::stringstream& ss, P p, Q...q)
{
    // Prepend 6 zeroes "0" to p and repeat for the other arguments. 
    ss << "-" << std::setfill('0') << std::setw(6) << p; 
    if constexpr (sizeof...(Q) > 0) 
        stream_labels(ss, q...); 
}

template<typename... Q> 
std::string vtk_file_name(const std::string& fileName, Q... q)
{
    std::stringstream ss; 
    ss << fileName; 
    stream_labels(ss, q...); 
    ss << ".vtk"; 
    return ss.str(); 
} 

template <typename Geometry, typename ... Q>
void write_vtk(
    Geometry const& geom, 
    std::string const& fileName, Q ... q
) // Usage example: write_vtk(geo, fileName, cellLabel, faceLabel, pointLabel); 
{
    // Append the labels q... paded with zeros to fileName.
    auto appendedName = vtk_file_name(fileName, q...); 
    write_vtk(geom, appendedName); 
}

} // End namespace geophase 

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //
