/********************************************************
 *  ██████╗  ██████╗████████╗██╗
 * ██╔════╝ ██╔════╝╚══██╔══╝██║
 * ██║  ███╗██║        ██║   ██║
 * ██║   ██║██║        ██║   ██║
 * ╚██████╔╝╚██████╗   ██║   ███████╗
 *  ╚═════╝  ╚═════╝   ╚═╝   ╚══════╝
 * Geophysical Computational Tools & Library (GCTL)
 *
 * Copyright (c) 2023  Yi Zhang (yizhang-geo@zju.edu.cn)
 *
 * GCTL is distributed under a dual licensing scheme. You can redistribute 
 * it and/or modify it under the terms of the GNU Lesser General Public 
 * License as published by the Free Software Foundation, either version 2 
 * of the License, or (at your option) any later version. You should have 
 * received a copy of the GNU Lesser General Public License along with this 
 * program. If not, see <http://www.gnu.org/licenses/>.
 * 
 * If the terms and conditions of the LGPL v.2. would prevent you from using 
 * the GCTL, please consider the option to obtain a commercial license for a 
 * fee. These licenses are offered by the GCTL's original author. As a rule, 
 * licenses are provided "as-is", unlimited in time for a one time fee. Please 
 * send corresponding requests to: yizhang-geo@zju.edu.cn. Please do not forget 
 * to include some description of your company and the realm of its activities. 
 * Also add information on how to contact you by electronic and paper mail.
 ******************************************************/

#include "meshdata.h"
#include "fstream"

gctl::meshdata::meshdata()
{
    name_ = "untitled";
    loctype_ = NodeData;
    valtype_ = Scalar;
    output_ok_ = true;
    nan_val_ = GCTL_BDL_MAX;
}

gctl::meshdata::meshdata(mesh_data_type_e in_loctype, 
    mesh_data_value_e in_valtype, size_t size, 
    std::string name, bool if_output, double nan_val)
{
    create(in_loctype, in_valtype, size, name, if_output, nan_val);
}

gctl::meshdata::~meshdata()
{
    clear();
}

void gctl::meshdata::clear()
{
    name_ = "untitled";
    loctype_ = NodeData;
    valtype_ = Scalar;
    output_ok_ = true;
    nan_val_ = GCTL_BDL_MAX;
    datval_.clear();
}

void gctl::meshdata::create(mesh_data_type_e in_loctype, 
    mesh_data_value_e in_valtype, size_t size,
    std::string name, bool if_output, double nan_val)
{
    name_ = name;
    loctype_ = in_loctype;
    valtype_ = in_valtype;
    output_ok_ = if_output;
    nan_val_ = nan_val;
    if (size != 0) datval_.resize(size, 0.0);
    return;
}

gctl::array<gctl::point3dc> gctl::meshdata::export_vector() const
{
    if (valtype_ != Vector) throw std::runtime_error("[gctl::meshdata::export_vector] Invalid data type.");

    size_t n = 0;
    array<point3dc> tmp(datval_.size()/3);
    for (size_t i = 0; i < datval_.size(); i += 3)
    {
        tmp[n] = point3dc(datval_[3*i], datval_[3*i + 1], datval_[3*i + 2]);
        n++;
    }
    return tmp;
}

gctl::array<gctl::tensor> gctl::meshdata::export_tensor() const
{
    if (valtype_ != Tensor) throw std::runtime_error("[gctl::meshdata::export_tensor] Invalid data type.");

    size_t n = 0;
    array<tensor> tmp(datval_.size()/9);
    for (size_t i = 0; i < datval_.size(); i += 9)
    {
        tmp[n] = tensor(datval_[9*i], datval_[9*i + 1], datval_[9*i + 2],
            datval_[9*i + 3], datval_[9*i + 4], datval_[9*i + 5],
            datval_[9*i + 6], datval_[9*i + 7], datval_[9*i + 8]);
        n++;
    }
    return tmp;
}

void gctl::meshdata::show_info(std::ostream &os) const
{
    os << "Data: " << name_ << " | ";
	if (loctype_ == NodeData) os << "node data | ";
	if (loctype_ == ElemData) os << "element data | ";
	if (loctype_ == ElemData2D) os << "2D element data | ";
	if (loctype_ == ElemData3D) os << "3D element data | ";
	if (valtype_ == Scalar) os << "scalar | ";
	if (valtype_ == Vector) os << "vector | ";
	if (valtype_ == Tensor) os << "tensor | ";
	if (output_ok_) os << "output" << std::endl;
	else os<< "no-output" << std::endl;
	return;
}

void gctl::meshdata::show_stats(std::ostream &os) const
{
    int vn = 0;
    if (valtype_ == Vector)
    {
        for (size_t i = 0; i < datval_.size(); i += 3)
        {
            if (fabs(datval_[3*i] - nan_val_) > 1e-10) vn++;
        }

        array<double> tmp_x(vn), tmp_y(vn), tmp_z(vn);

        vn = 0;
        for (size_t i = 0; i < datval_.size(); i += 3)
        {
            if (fabs(datval_[3*i] - nan_val_) > 1e-10)
            {
                tmp_x[vn] = datval_[3*i];
                tmp_y[vn] = datval_[3*i + 1];
                tmp_z[vn] = datval_[3*i + 2];
                vn++;
            }
        }
        
        os  << "max  = (" << tmp_x.max() << ", " << tmp_y.max() << ", " << tmp_z.max() << ")\n"
            << "min  = (" << tmp_x.min() << ", " << tmp_y.min() << ", " << tmp_z.min() << ")\n"
            << "mean = (" << tmp_x.mean() << ", " << tmp_y.mean() << ", " << tmp_z.mean() << ")\n"
            << "std  = (" << tmp_x.std() << ", " << tmp_y.std() << ", " << tmp_z.std() << ")\n"
            << "rms  = (" << tmp_x.rms() << ", " << tmp_y.rms() << ", " << tmp_z.rms() << ")\n";
    }
    else if (valtype_ == Tensor)
    {
        for (size_t i = 0; i < datval_.size(); i += 9)
        {
            if (fabs(datval_[9*i] - nan_val_) > 1e-10) vn++;
        }

        array<double> tmp_xx(vn), tmp_xy(vn), tmp_xz(vn);
        array<double> tmp_yx(vn), tmp_yy(vn), tmp_yz(vn);
        array<double> tmp_zx(vn), tmp_zy(vn), tmp_zz(vn);

        vn = 0;
        for (size_t i = 0; i < datval_.size(); i += 9)
        {
            if (fabs(datval_[9*i] - nan_val_) > 1e-10)
            {
                tmp_xx[vn] = datval_[9*i];
                tmp_xy[vn] = datval_[9*i + 1];
                tmp_xz[vn] = datval_[9*i + 2];
                tmp_yx[vn] = datval_[9*i + 3];
                tmp_yy[vn] = datval_[9*i + 4];
                tmp_yz[vn] = datval_[9*i + 5];
                tmp_zx[vn] = datval_[9*i + 6];
                tmp_zy[vn] = datval_[9*i + 7];
                tmp_zz[vn] = datval_[9*i + 8];
                vn++;
            }
        }

        os << "max = \n"
        << tmp_xx.max() << ", " << tmp_xy.max() << ", " << tmp_xz.max() << "\n"
        << tmp_yx.max() << ", " << tmp_yy.max() << ", " << tmp_yz.max() << "\n"
        << tmp_zx.max() << ", " << tmp_zy.max() << ", " << tmp_zz.max() << "\n";

        os << "min = \n"
        << tmp_xx.min() << ", " << tmp_xy.min() << ", " << tmp_xz.min() << "\n"
        << tmp_yx.min() << ", " << tmp_yy.min() << ", " << tmp_yz.min() << "\n"
        << tmp_zx.min() << ", " << tmp_zy.min() << ", " << tmp_zz.min() << "\n";

        os << "mean = \n" 
        << tmp_xx.mean() << ", " << tmp_xy.mean() << ", " << tmp_xz.mean() << "\n"
        << tmp_yx.mean() << ", " << tmp_yy.mean() << ", " << tmp_yz.mean() << "\n"
        << tmp_zx.mean() << ", " << tmp_zy.mean() << ", " << tmp_zz.mean() << "\n";

        os << "std  = \n"
        << tmp_xx.std() << ", " << tmp_xy.std() << ", " << tmp_xz.std() << "\n"
        << tmp_yx.std() << ", " << tmp_yy.std() << ", " << tmp_yz.std() << "\n"
        << tmp_zx.std() << ", " << tmp_zy.std() << ", " << tmp_zz.std() << "\n";
        
        os << "rms  = \n"
        << tmp_xx.rms() << ", " << tmp_xy.rms() << ", " << tmp_xz.rms() << "\n"
        << tmp_yx.rms() << ", " << tmp_yy.rms() << ", " << tmp_yz.rms() << "\n"
        << tmp_zx.rms() << ", " << tmp_zy.rms() << ", " << tmp_zz.rms() << "\n";
    }
    else // Scalar
    {
        for (size_t i = 0; i < datval_.size(); i++)
        {
            if (fabs(datval_[i] - nan_val_) > 1e-10) vn++;
        }

        array<double> tmp(vn);

        vn = 0;
        for (size_t i = 0; i < datval_.size(); i++)
        {
            if (fabs(datval_[i] - nan_val_) > 1e-10)
            {
                tmp[vn] = datval_[i];
                vn++;
            }
        }
        
        os << "max = " << tmp.max() << ", min = " << tmp.min() << std::endl;
        os << "mean = " << tmp.mean() << ", std = " << tmp.std() << ", rms = " << tmp.rms() << std::endl;
    }
	return;
}

void gctl::meshdata::load_binary(std::ifstream &infile)
{
    int name_size;
    infile.read((char*) &name_size, sizeof(int));

    name_.resize(name_size);
    infile.read((char*) name_.c_str(), name_size);

    infile.read((char*) &loctype_, sizeof(int));
	infile.read((char*) &valtype_, sizeof(int));
    infile.read((char*) &output_ok_, sizeof(bool));

    int n;
    infile.read((char*) &n, sizeof(int));
    datval_.resize(n);

	infile.read((char*) datval_.get(), sizeof(double)*datval_.size());
    infile.read((char*) &nan_val_, sizeof(double));
    return;
}

void gctl::meshdata::save_binary(std::ofstream &outfile)
{
    int name_size = name_.size();
	outfile.write((char*) &name_size, sizeof(int));
	outfile.write((char*) name_.c_str(), name_size);
	outfile.write((char*) &loctype_, sizeof(int));
	outfile.write((char*) &valtype_, sizeof(int));
	outfile.write((char*) &output_ok_, sizeof(bool));

    int n = datval_.size();
	outfile.write((char*) &n, sizeof(int));
	outfile.write((char*) datval_.get(), sizeof(double)*datval_.size());
    outfile.write((char*) &nan_val_, sizeof(double));
    return;
}