#include "../src/nmie-basic.hpp"
#include "../src/nmie-nearfield.hpp"
#include "gtest/gtest.h"
#include "test_cases.hpp"

// TEST(RunFieldCalculationCartesian, DISABLED_HandlesInput) {
TEST(RunFieldCalculationCartesian, HandlesInput) {
  nmie::MultiLayerMie<nmie::FloatType> nmie;
  //  EXPECT_THROW(nmie.RunFieldCalculationPolar(0), std::invalid_argument);
  //  EXPECT_THROW(nmie.RunFieldCalculationPolar(1,1,10,5),
  //  std::invalid_argument);
  nmie::FloatType total_r = 2 * nmie::PI_ * 1000 / 532;
  //  double r = 1500;
  nmie.SetLayersSize({total_r / 2, total_r});
  nmie.SetLayersIndex({{1.330, 0}, {1.33, 0}});
  nmie.RunMieCalculation();
  double relative_max_distance = 1e-10;
  //  nmie.SetModeNmaxAndType(3,-1);
  //  int nmax = 21;
  nmie.RunFieldCalculationCartesian(2, 5, relative_max_distance,
                                    nmie::Planes::kEk, 1.0, 0, 0, false, 3);
  auto Eabs = nmie.GetFieldEabs();
  auto E = nmie.GetFieldE();
  std::cout << std::endl;
  {
    // Eabs points are located near the sphere outer border
    //
    //    0   1   2   3   4
    //    ----- border ----
    //    5   6   7   8   9
    // distance between points (0) and (4) is relative_max_distance*total_r,
    // initial value used for the test was 1e-10*total_r, so we expect good
    // linear dependence for points from 0 to 4 and 5 to 9. In the asserts we
    // check, that the slope doesn't change too fast inside the curve. While
    // writing this, the test was failing. The value of z-coordinates of 2 and 7
    // points = 0
    using nmie::nmm::abs;
    EXPECT_TRUE((abs(Eabs[0] - Eabs[1]) + abs(Eabs[3] - Eabs[4])) >=
                abs(Eabs[1] - Eabs[2]));
    EXPECT_TRUE((abs(Eabs[5] - Eabs[6]) + abs(Eabs[8] - Eabs[9])) >=
                abs(Eabs[6] - Eabs[7]));
  }

  //  nmie.RunFieldCalculationCartesian(2, 2, 2, nmie::Planes::kHk,
  //                                    0, 0, 0, true);
  //  nmie.RunFieldCalculationCartesian(2, 2, 2, nmie::Planes::kEH,
  //                                    0, 0, 0, true);
  // TODO add check of E and H symmetry for X and Y axis inversion

  //  EXPECT_EQ(1, nmie.GetMaxTerms());
  //  EXPECT_FALSE(nmie.GetFieldConvergence());
  //  auto Eabs = nmie.GetFieldEabs();
  //  EXPECT_TRUE(std::isnan(static_cast<double>(Eabs[0])));
}
TEST(LargeBubbleSpectrum, DISABLED_HandlesInput) {
  // TEST(LargeBubbleSpectrum, HandlesInput) { // TODO fix fail...
  nmie::MultiLayerMie<nmie::FloatType> nmie;
  nmie::FloatType core_r = 2 * nmie::PI_ * 100;
  nmie::FloatType shell_r = 2 * nmie::PI_ * (100 + 0.1);
  nmie.SetLayersIndex({{1, 0}, {1.33, 0}});
  double central_WL = 0.7007;
  double relative_distance = 1e-10;
  double dWL = central_WL * relative_distance;
  std::vector<double> Qsca(5);
  for (int i = 0; i < 5; ++i) {
    auto WL = static_cast<nmie::FloatType>(central_WL + (i - 2) * dWL);
    nmie.SetLayersSize({core_r / WL, shell_r / WL});
    nmie.RunMieCalculation();
    Qsca[i] = static_cast<double>(nmie.GetQsca());
    std::cout << "Qsca[" << i << "]=" << Qsca[i] << std::endl;
  }

  {
    // Eabs points are located near the sphere outer border
    //
    //    0   1   2   3   4
    //    ------- WL ------>
    // distance between points (0) and (4) is 5*relative_distance*central_WL,
    // initial value used for the test was 5*1e-10*0.7007, so we expect good
    // linear dependence for points from 0 to 4. In the asserts we check, that
    // the slope doesn't change too fast inside the curve.
    using std::abs;
    EXPECT_TRUE((abs(Qsca[0] - Qsca[1]) + abs(Qsca[3] - Qsca[4])) >=
                abs(Qsca[1] - Qsca[2]));
  }
}

// TEST(RunFieldCalculationPolar, DISABLED_HandlesInput) {
TEST(RunFieldCalculationPolar, HandlesInput) {
  nmie::MultiLayerMie<nmie::FloatType> nmie;
  EXPECT_THROW(nmie.RunFieldCalculationPolar(0), std::invalid_argument);
  EXPECT_THROW(nmie.RunFieldCalculationPolar(1, 1, 10, 5),
               std::invalid_argument);
  double r = 60;
  //  double r = 1500;
  nmie.SetLayersSize({r / 2, r});
  nmie.SetLayersIndex({{1.33, 0}, {1.33, 0}});
  nmie.RunMieCalculation();
  nmie.RunFieldCalculationPolar(1, 1, 0.5145 * r, r * 0.5148, 0, 3.14, 0, 0,
                                true, 1);
  EXPECT_EQ(1, nmie.GetMaxTerms());
  EXPECT_FALSE(nmie.GetFieldConvergence());
  auto Eabs = nmie.GetFieldEabs();
  EXPECT_TRUE(std::isnan(static_cast<double>(Eabs[0])));
}
//#ifndef MULTI_PRECISION
// TEST(BulkSphere, DISABLED_HandlesInput) {
TEST(BulkSphere, HandlesInput) {
  nmie::MultiLayerMie<nmie::FloatType> nmie;
  for (const auto& data : parameters_bulk_sphere) {
    auto x = std::get<0>(data);
    auto m = std::get<1>(data);
    nmie.SetLayersSize({x});
    nmie.SetLayersIndex({m});
    nmie.SetMaxTerms(-1);
    //    nmie.RunMieCalculation();
    //    std::cout<<" test case: "<<std::get<2>(data)<<"
    //    Qsca="<<nmie.GetQsca()<<std::endl;
    nmie.RunFieldCalculationPolar(4, 3, x, x * 3, 0,
                                  static_cast<double>(nmie::PI_), 0,
                                  static_cast<double>(nmie::PI_), true, -1);
    auto Eabs = nmie.GetFieldEabs();
    for (auto& E : Eabs)
      E = nmie::pow2(E);
    //    print(Eabs)
    EXPECT_TRUE(nmie.GetFieldConvergence())
        << "Outside test for x=" << x << " m=" << m
        << " test case: " << std::get<2>(data) << std::endl;
    nmie.RunFieldCalculationPolar(4, 10, x * 0.01, x, 0,
                                  static_cast<double>(nmie::PI_), 0,
                                  static_cast<double>(nmie::PI_), true, -1);
    EXPECT_TRUE(nmie.GetFieldConvergence())
        << "Inside test for x=" << x << " m=" << m
        << " test case: " << std::get<2>(data) << std::endl;
  }
}
//#endif

int main(int argc, char** argv) {
  testing::InitGoogleTest(&argc, argv);
  return RUN_ALL_TESTS();
}