/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | faSavageHutterFOAM
    \\  /    A nd           | Copyright (C) 2017 Matthias Rauter
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM 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.

    OpenFOAM 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 OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Application
    fafiedlsToAscii

Description
    Write finite area fields into a csv file to create raster with gdal

Author
    Matthias Rauter matthias.rauter@uibk.ac.at

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


/*
 * convert to raster:
 *
 *  add file faFields.vrt containing:
 * <OGRVRTDataSource>
 *    <OGRVRTLayer name="faFields">
 *        <SrcDataSource>faFields.csv</SrcDataSource>
 *        <GeometryType>wkbPoint</GeometryType>
 *        <LayerSRS>WGS84</LayerSRS>
 *        <GeometryField encoding="PointFromColumns" x="x" y="y"/>
 *    </OGRVRTLayer>
 * </OGRVRTDataSource>
 *
 * and run:
 *
 * gdal_grid -a nearest:radius1=10:radius2=10:nodata=-9999 -outsize 1000 2000 -zfield z -l faFields faFields.vrt z.tiff
 *
 *
 *
 */

#include "fvCFD.H"
#include "faCFD.H"
#include "IOobjectList.H"

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

int main(int argc, char *argv[])
{
#   include "setRootCase.H"
#   include "createTime.H"
#   include "createMesh.H"
#   include "createFaMesh.H"

    areaVectorField c = aMesh.areaCentres();

    IOdictionary fafieldsToAsciiDict
    (
         IOobject
         (
              "fafieldsToAsciiDict",
              runTime.system(),
              mesh,
              IOobject::MUST_READ,
              IOobject::NO_WRITE
         )
    );


    Foam::instantList timeDirs = Foam::timeSelector::select0(runTime, args);

    List<word> fields;
    List<areaVectorField*> vfl;
    List<areaScalarField*> sfl;

    vector offset;
    fafieldsToAsciiDict.lookup("offset") >> offset;

    fafieldsToAsciiDict.lookup("fields") >> fields;


    forAll(timeDirs, timeI)
    {
        forAll(vfl, i)
            delete vfl[i];
        vfl.resize(0);
        forAll(sfl, i)
            delete sfl[i];
        sfl.resize(0);

        runTime.setTime(timeDirs[timeI], timeI);
        Foam::Info << "\nTime = " << runTime.timeName() << endl;

        mesh.readUpdate();

        IOobjectList objects(mesh, runTime.timeName());


        std::ostringstream buf;
        buf << runTime.timeName();
        buf << "/faFields.csv";

        std::ofstream csvfile;
        csvfile.open(buf.str().c_str());
        csvfile << "x;" << "y;" << "z;";

        forAll(fields, fi)
        {
            if (objects.lookup(fields[fi])->headerClassName() == "areaVectorField")
            {
                Info << "reading areaVectorField " << objects.lookup(fields[fi])->name() << endl;

                areaVectorField *vf =
                    new areaVectorField(
                    IOobject
                    (
                        fields[fi],
                        runTime.timeName(),
                        mesh,
                        IOobject::MUST_READ,
                        IOobject::NO_WRITE
                    ),
                    aMesh
                );
                vfl.resize(vfl.size()+1);
                vfl[vfl.size()-1] = vf;

                csvfile << fields[fi] << "x;" << fields[fi] << "y;" << fields[fi] << "z" ;
            }
            else if (objects.lookup(fields[fi])->headerClassName() == "areaScalarField")
            {
                Info << "reading areaScalarField " << objects.lookup(fields[fi])->name() << endl;

                areaScalarField *sf =
                    new areaScalarField
                    (
                        IOobject
                        (
                            fields[fi],
                            runTime.timeName(),
                            mesh,
                            IOobject::MUST_READ,
                            IOobject::NO_WRITE
                        ),
                        aMesh
                    );
                sfl.resize(sfl.size()+1);
                sfl[vfl.size()-1] = sf;

                csvfile << fields[fi];
            }
            else
            {
                Info << "object" << objects.lookup(fields[fi])->name() << "is neighter areaScalarField nor areaVectorField!" << endl;
            }
            if (fi != fields.size()-1)
            {
                csvfile << ";";
            }
            else
            {
                csvfile << nl;
            }
        }
        forAll(aMesh.areaCentres(), i)
        {
            csvfile << aMesh.areaCentres()[i][0]-offset[0] << ";" << aMesh.areaCentres()[i][1]-offset[1] << ";" << aMesh.areaCentres()[i][2]-offset[2] << ";";

            label sfi = 0;
            label vfi = 0;
            forAll(fields, fi)
            {
                if (objects.lookup(fields[fi])->headerClassName() == "areaVectorField")
                {
                    areaVectorField *vf = vfl[vfi];
                    csvfile << (*vf)[i][0] << ";" << (*vf)[i][1] << ";" << (*vf)[i][2];
                    vfi++;
                }
                else if (objects.lookup(fields[fi])->headerClassName() == "areaScalarField")
                {
                    areaScalarField *sf = sfl[sfi];
                    csvfile << (*sf)[i];
                    sfi++;
                }
                if (fi != fields.size()-1)
                {
                    csvfile << ";";
                }
                else
                {
                    csvfile << nl;
                }
            }

        }


        csvfile.close();

    }

    return(0);

}

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