testTabulatedFunction.cc

Go to the documentation of this file.

#include <utl/TabulatedFunction.h>
#include <utl/TabulatedFunctionErrors.h>
#include <utl/MultiTabulatedFunction.h>
#include <utl/MultiTabulatedFunctionErrors.h>
#include <utl/AugerException.h>
#include <iostream>
#include <iomanip>
#include <cstdlib>
#include <fstream>

using namespace std;

using namespace utl;


int
main()
{
  const double aX[] = { 1.1, 2.3, 4.2, 4.5,  6.7, 8.333, 12.12, 13.45 };
  const double aY[] = { 5.7, 3.3, 6.3, 11.1, 1.1, -18.4, -20.7, 44.3 };
  const int aN = 8;

  const double a2X[] = { 1, 2, 4, 3, 5, 6, 7, 8, 9, 10 };
  const double a2Y[] = { 1, 2, 4, 3, 5, 6, 7, 8, 9, 10 };
  const int a2N = 10;

  const std::vector<double> vecX(aX, aX + aN);
  const std::vector<double> vecY(aY, aY + aN);
  const std::vector<double> vec2X(a2X, a2X + a2N);
  const std::vector<double> vec2Y(a2Y, a2Y + a2N);

  TabulatedFunction TF(vecX, vecY);
  TabulatedFunction TF2(vec2X, vec2Y);

  const double tx[] = { 2.11, 5.6, 2.7, 8.333, 9 };
  const int tN = 5;
  cout << "Test points TF:\n";
  for (int i = 0; i < tN; ++i)
    cout << "x= " << setw(5) << tx[i] << "  y= " << setw(5) << TF.Y(tx[i]) << '\n';
  cout << endl;

  const double t2x[] = { 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5 };
  const int t2N = 7;
  cout << "Test points TF2:\n";
  for (int i = 0; i < t2N; ++i)
    cout << "x= " << setw(3) << t2x[i] << "  y= " << setw(5) << TF2.Y(t2x[i]) << '\n';
  cout << endl;

  double a3X[] = { 1, 3, 5 };
  double a3Y[] = { 1, 9, 25 };
  const int a3N = 3;
  const std::vector<double> vec3X(a3X, a3X + a3N);
  const std::vector<double> vec3Y(a3Y, a3Y + a3N);
  const TabulatedFunction TF3(vec3X, vec3Y);

  cout << "Test points TF3:\n";
  for (double xVal = 1.0; xVal <= 5.0; xVal += 0.1)
    cout << "x= " << setw(3) << xVal << "  y= " << setw(5) << TF3.Y(xVal) << '\n';
  cout << endl;

  cout << "Test points TF3 - interpolation order 1:\n";
  for (double xVal = 1.0; xVal <= 5.0 ; xVal += 0.1)
    cout << "x= " << setw(3) << xVal << "  y= " << setw(5) << TF3.InterpolateY(xVal, 1) << '\n';
  cout << endl;

  cout << "Test points TF3 - interpolation order 2:\n";
  for (double xVal = 1.0; xVal <= 5.0 ; xVal += 0.1)
    cout << "x= " << setw(3) << xVal << "  y= " << setw(5) << TF3.InterpolateY(xVal, 2) << '\n';

  cout << endl << "Testing operator[]..." << endl
       << TF[0].X() << " " << TF[0].Y() << endl
       << endl << "Testing iterator" << endl;

  utl::TabulatedFunction::Iterator it;

  for (utl::TabulatedFunction::Iterator it = TF.Begin();
       it != TF.End(); ++it)
    cout << it->X() << " " << it->Y() << endl;

  cout << endl << "Testing operator[]..." << endl;
  for (unsigned int i = 0; i < TF.GetNPoints(); ++i)
    cout << TF[i].X() << " " << TF[i].Y() << endl;

  cout << endl
       << endl << "Testing PushBack..." << endl;

  TabulatedFunction f;
  f.PushBack(1, 2);    
  f.PushBack(3, 4);  
  f.PushBack(5, 6);
  f.PushBack(7, 8);   
  f.PushBack(9, 10);

  for (utl::TabulatedFunction::Iterator it = f.Begin();
       it != f.End(); ++it)
    cout << it->X() << " " << it->Y() << endl;

  cout << "Bounds " << f.GetX(0) << " " << f.GetX(f.GetNPoints()-1) << endl;
  cout << endl << "Testing const iterator" << endl;

  const TabulatedFunction cf(vecX, vecY);

  for (utl::TabulatedFunction::ConstIterator it = cf.Begin();
       it != cf.End(); ++it)
    cout << it->X() << " " << it->Y() << endl;
  cout << endl;

  utl::TabulatedFunction tab;
  {   
    const double a = 1;
    double b = 2;
    double c = 3;
    const double d = 4;

    tab.PushBack(a, b);
    tab.PushBack(c, d);
  }

  for (unsigned int i = 0; i < tab.GetNPoints(); ++i)
    cout << tab[i].X() << " " << tab[i].Y() << endl;

  cout << "Testing TabulatedFunctionError" << endl;

  TabulatedFunctionErrors tberr;

  tberr.PushBack(11., 0.5, 1., 0.1);
  tberr.PushBack(12., 0.6, 2., 0.2);
  tberr.PushBack(13., 0.7, 3., 0.3);
  tberr.PushBack(14., 0.8, 4., 0.4);    
  tberr.PushBack(15., 0.9, 5., 0.5);
  tberr.PushBack(16., 1.0, 6., 0.6);
  tberr.PushBack(17., 1.1, 7., 0.7);

  cout << " iterator " << endl;

  for (TabulatedFunctionErrors::IteratorErr it = tberr.Begin();
       it != tberr.End(); ++it)
    cout << "x " << it->X() << " +- " << it->XErr()
         << ", y " << it->Y() << " +- " << it->YErr() << endl;
  cout << endl;

  cout << "Y(16.5)= " << tberr.Y(16.5) << endl;
  cout << "Y(11.5)= " << tberr.Y(11.5) << endl;

  for (int i = 0; i < 3; ++i)
     cout << "x " << tberr[i].X() << " +- " << tberr[i].XErr()
          << ", y " << tberr[i].Y() << " +- " << tberr[i].YErr() << endl;

  cout << "Test MultiTabulatedFunction..." << endl;

  MultiTabulatedFunction mtabf;
  MultiTabulatedFunctionErrors mtabferr;

  mtabf.AddTabulatedFunction(f, 1);
  mtabferr.AddTabulatedFunctionErrors(tberr, 2);

  TabulatedFunction& fp = mtabf.GetTabulatedFunction(1);
  cout << "x " << fp[0].X() <<endl;

  TabulatedFunctionErrors& fp2 = mtabferr.GetTabulatedFunctionErrors(2);
  cout << "x " << fp2[0].X() << " +- " << fp2[0].XErr() << endl;

  cout << "Test copy ctor" << endl;
  MultiTabulatedFunction copy = mtabf;
  TabulatedFunction& fpcopy1 = copy.GetTabulatedFunction(1);
  cout << "x " << fpcopy1[0].X() << endl;
   
  MultiTabulatedFunctionErrors copye = mtabferr;
  TabulatedFunctionErrors& fpcopy2 = copye.GetTabulatedFunctionErrors(2);
  cout << "x " << fpcopy2[0].X() << " +- " << fpcopy2[0].XErr() << endl;

  // see what happens if you initialize a tabluated function with just one pair
  const std::vector<double> sX(1, 1.0);
  const std::vector<double> sY(1, 10.0);
  TabulatedFunction STF(sX, sY);
  cout << "tabulated func initialzed with one pair" << endl;

  try {
    cout << " STF(1.0) " << STF.Y(1.0) << endl;
  } catch (AugerException& ex) {
    cout << "message " << ex.GetMessage() << endl;
  }

  try {
    cout << " STF(0.5) " << STF.Y(0.5) << endl;
  } catch (AugerException& ex) {
    cout << "caught expected exception " << ex.GetExceptionName() << endl;
  }

  cout << "checking FindX and FindY methods" << endl;
  for(TabulatedFunctionErrIterator it = tberr.Begin(); it != tberr.End(); it++)
  {
      cout << it->X() << " == " << tberr.FindY(it->Y())->X() << "?\t" 
           << it->Y() << " == " << tberr.FindX(it->X())->Y() << "?" << endl;
  }

  // try to push back two identical points
  //
  cout << "pushing the same point twice " << endl;
  TabulatedFunction crap;
  crap.PushBack(5.,  15.);
  crap.PushBack(10., 20.);
  crap.PushBack(10., 20.);

  cout << "pushed 3 points, got " << crap.GetNPoints() << '\n'
       << "crap(7) = " << crap.InterpolateY(7., 3) << endl;

  return EXIT_SUCCESS;
}


// Configure (x)emacs for this file ...
// Local Variables:
// mode: c++
// compile-command: "make -C .. testTabulatedFunction -k"
// End: