FloatFun.C

Go to the documentation of this file.

#include "jbltools/sfh/FloatFun.h"
#include "jbltools/sfh/FloatFunHelpers.h"
#include "jbltools/sfh/IntFun.h"

#include "jbltools/sfh/FillIterator.h"

#include <cmath>
// #include <limits>
#include <values.h>
#include <iostream>
#include <cassert>

#include <stdio.h>
#include <cstring>


FloatFun& abs (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::abs>("abs", arg_);
}

FloatFun& acos (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::acos>("acos", arg_);
}

// FloatFun& acosh (const FloatFun& arg_) {
//   return *new AnyFloatFun1<std::acosh>("acosh", arg_);
// }
// 
FloatFun& asin (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::asin>("asin", arg_);
}

// FloatFun& asinh (const FloatFun& arg_) {
//   return *new AnyFloatFun1<std::asinh>("asinh", arg_);
// }
// 
FloatFun& atan (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::atan>("atan", arg_);
}

// FloatFun& atanh (const FloatFun& arg_) {
//   return *new AnyFloatFun1<std::atanh>("atanh", arg_);
// }
// 
FloatFun& ceil (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::ceil>("ceil", arg_);
}

FloatFun& cos (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::cos>("cos", arg_);
}

FloatFun& cosh (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::cosh>("cosh", arg_);
}

FloatFun& exp (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::exp>("exp", arg_);
}

FloatFun& fabs (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::fabs>("fabs", arg_);
}

FloatFun& floor (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::floor>("floor", arg_);
}

FloatFun& log (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::log>("log", arg_);
}

FloatFun& log10 (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::log10>("log10", arg_);
}

FloatFun& sin (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::sin>("sin", arg_);
}

FloatFun& sinh (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::sinh>("sinh", arg_);
}

FloatFun& sqrt (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::sqrt>("sqrt", arg_);
}

FloatFun& tan (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::tan>("tan", arg_);
}

FloatFun& tanh (const FloatFun& arg_) {
  return *new AnyFloatFun1<std::tanh>("tanh", arg_);
}

FloatFun& square (const FloatFun& arg_) {
  return *new AnyFloatFun1<square>("square", arg_);
}

FloatFun& rad (const FloatFun& arg_) {
  return *new AnyFloatFun1<rad>("rad", arg_);
}

FloatFun& deg (const FloatFun& arg_) {
  return *new AnyFloatFun1<deg>("deg", arg_);
}

FloatFun& atan2 (const FloatFun& lhs, const FloatFun& rhs) {
  return *new AnyFloatFun2<std::atan2>("atan2", lhs, rhs);
}

FloatFun& pow (const FloatFun& arg1, const FloatFun& arg2) {
  return *new AnyFloatFun2<std::pow>("pow", arg1, arg2);
}

Float_FF MaxFloatFun::operator() () const {
  // Float_FF result = std::numeric_limits<Float_FF>::min();
  Float_FF result = -(MAXFLOAT);
  if (intiter->reset()) {
    do {
      if (!cut || (*cut)()) {
        Float_FF f = fun ? (*fun)() : 0;
        if (f > result) result = f;
      }
    } while (intiter->next());
  }
  return result;
}             

Float_FF MinFloatFun::operator() () const {
//  Float_FF result = std::numeric_limits<Float_FF>::max();
  Float_FF result = MAXFLOAT;
  if (intiter->reset()) {
    do {
      if (!cut || (*cut)()) {
        Float_FF f = fun ? (*fun)() : 0;
        if (f < result) result = f;
      }
    } while (intiter->next());
  }
  return result;
}             

Float_FF SumFloatFun::operator() () const {
  Float_FF sum = 0;
  if (intiter->reset()) {
    do {
      if (!cut || (*cut)()) {
        Float_FF f = fun ? (*fun)() : 1;
        sum += f;
      }
    } while (intiter->next());
  }
  return sum;
}             

Float_FF AverageFloatFun::operator() () const {
  Float_FF sum = 0;
  int   n = 0;
  if (intiter->reset()) {
    do {
      if (!cut || (*cut)()) {
        Float_FF f = fun ? (*fun)() : 1;
        sum +=f;
        ++n;
      }
    } while (intiter->next());
  }
  return n ? (sum/n) : 0;
}             

OneFunCompound::OneFunCompound (const std::string& funname_, const FloatFun& arg_)
: FloatFun (funname_+'('+arg_.getName()+')'), arg (arg_)
{}

OneFunCompound::OneFunCompound (const char *funname_, const FloatFun& arg_)
: FloatFun (std::string (funname_)+'('+arg_.getName()+')'), arg (arg_)
{}

OneFunCompound::OneFunCompound (const FloatFun& arg_, const std::string& name_)
: FloatFun (name_), arg (arg_)
{}
   
TwoFunCompound::TwoFunCompound (const FloatFun& lhs_, const std::string opname_, const FloatFun& rhs_) 
: FloatFun (str('(')+lhs_.getName()+')'+opname_+'('+rhs_.getName()+')'), lhs (lhs_), rhs (rhs_)
{}

TwoFunCompound::TwoFunCompound (const std::string funname_, const FloatFun& lhs_, const FloatFun& rhs_) 
: FloatFun (funname_+'('+lhs_.getName()+','+rhs_.getName()+')'), lhs (lhs_), rhs (rhs_)
{}

TwoFunCompound::TwoFunCompound (const FloatFun& lhs_, const FloatFun& rhs_, const std::string name_) 
: FloatFun (name_), lhs (lhs_), rhs (rhs_)
{}

const FillIterator *TwoFunCompound::getIterator() const { 
  assert (!lhs.getIterator() || !rhs.getIterator() || lhs.getIterator()==rhs.getIterator());
  return (lhs.getIterator()) ? lhs.getIterator() : rhs.getIterator();
}

TwoFloatFunComp::TwoFloatFunComp (const FloatFun& lhs_, const std::string& opname_, const FloatFun& rhs_)
: FloatFunComp(str('(')+lhs_.getName()+')'+opname_+'('+rhs_.getName()+')'),
  lhs (lhs_), 
  rhs (rhs_) 
{
  assert (&lhs);
  assert (&rhs);
}
   
const FillIterator *TwoFloatFunComp::getIterator() const { 
  assert (!lhs.getIterator() || !rhs.getIterator() || lhs.getIterator()==rhs.getIterator());
  return (lhs.getIterator()) ? lhs.getIterator() : rhs.getIterator();
}

ChainComparison::ChainComparison (const FloatFunComp& lhs_, const std::string& opname_, const FloatFun& rhs_)
: FloatFunComp (lhs_.getName()+opname_+'('+rhs_.getName()+')'),
  lhs (lhs_), 
  rhs (rhs_)
{
  assert (&lhs);
  assert (&rhs);
}
    
const FillIterator *ChainComparison::getIterator() const { 
  assert (!lhs.getIterator() || !rhs.getIterator() || lhs.getIterator()==rhs.getIterator());
  return (lhs.getIterator()) ? lhs.getIterator() : rhs.getIterator();
}


ChainCompFloat::ChainCompFloat (const FloatFunComp& lhs_, const std::string& opname_, Float_FF rhs_)
: FloatFunComp (lhs_.getName()+opname_+str(rhs_)),
  lhs (lhs_), 
  rhs (rhs_)
{
  assert (&lhs);
}

FloatFloatFunComp::FloatFloatFunComp (Float_FF lhs_, const std::string& opname_, const FloatFun& rhs_)
: FloatFunComp (str(lhs_)+opname_+'('+rhs_.getName()+')'),
  lhs (lhs_), 
  rhs (rhs_) 
{
  assert (&rhs);
}

FloatFunFloatComp::FloatFunFloatComp (const FloatFun& lhs_, const std::string&  opname_, Float_FF rhs_) 
: FloatFunComp (str('(')+lhs_.getName()+')'+opname_+str(rhs_)),
  lhs (lhs_), 
  rhs (rhs_)
{
  assert (&lhs);
}

CollectiveFloatFun::CollectiveFloatFun (const std::string&           funname_,    
                                        const FloatFunPoRConst&      fun_,     
                                        const BaseCutPoRConst&       cut_,    
                                        const FillIteratorPoR&       intiter_,
                                        const FillIteratorPoRConst&  extiter_)
: FloatFun (makeName (funname_, fun_, cut_, intiter_, extiter_)),
  fun (fun_.pff), 
  cut (cut_.pbc), 
  intiter (intiter_.pfi), 
  extiter (extiter_.pfi)
{ 
  checkIterators();
  setName (makeName (funname_, fun, cut, intiter, extiter));
}       
                
std::string CollectiveFloatFun::makeName(const std::string&           funname_,       
                                         const FloatFunPoRConst&      fun_,        
                                         const BaseCutPoRConst&       cut_,      
                                         const FillIteratorPoR&       intiter_,  
                                         const FillIteratorPoRConst&  extiter_   
                                        ){
  const char *fname       = fun_.pff     ? fun_.pff->getName() : 0;
  const char *cutname     = cut_.pbc     ? cut_.pbc->getName() : 0;
  const char *intitername = intiter_.pfi ? intiter_.pfi->getName() : 0;
  const char *extitername = extiter_.pfi ? extiter_.pfi->getName() : 0;
  if (!fname)       fname="1";
  if (!cutname)     cutname="true";
  std::string result = funname_+'('+fname+", "+cutname+')';
  if (extitername) result = result + '['+intitername+']';
  
  return result;
}
             
void CollectiveFloatFun::checkIterators() {
  // Rules:
  // If no internal iterator is given, try to get iterator from
  // fun or cut; if both have an iterator, they must be equal or 
  // one of them must be equal to the external ietartor
  
  FillIterator *funiter = fun ? const_cast<FillIterator *>(fun->getIterator()) : 0;
  FillIterator *cutiter = cut ? const_cast<FillIterator *>(cut->getIterator()) : 0;
  
  if (!intiter) {
    // intiter not explicitly given, try to get it
    if (funiter && cutiter) {
      // funiter and cutiter defined
      if (funiter == cutiter) {
        if (funiter != extiter) {
          intiter = funiter;
          #if (DEBUG>=2)
          std::cerr << "CollectiveFloatFun '" << getName() << "': taking intiter from funiter '"
                    << funiter->getName() << "' == cutiter\n";
          #endif
        }
      }
      else {
        // funiter and cutiter different: One  must be equal to extiter,
        // the other should become intiter.
        if (funiter == extiter) {
          intiter = cutiter;
          std::cerr << "CollectiveFloatFun '" << getName() << "': taking intiter from cutiter '"
                    << cutiter->getName() << "' != funiter\n";
        }
        else if  (cutiter == extiter) {
          intiter = funiter;
          std::cerr << "CollectiveFloatFun '" << getName() << "': taking intiter from funiter '"
                    << funiter->getName() << "'!= cutiter\n";
        }
      }
    }
    else if (funiter) {
      // funiter defined, cutiter == 0
      if (funiter != extiter) {
        std::cerr << "CollectiveFloatFun '" << getName() << "': taking intiter from funiter'"
                    << funiter->getName() << "', cutiter==0\n";
        intiter = funiter;
      }
    }
    else if (cutiter) {
      // cutiter defined, funiter == 0
      if (cutiter != extiter) {
        std::cerr << "CollectiveFloatFun '" << getName() << "': taking intiter from cutiter'"
                    << cutiter->getName() << "', cutiter==0\n";
        intiter = cutiter;
      }
    }
    else {
      // funcutiter and cutiter == 0
      std::cerr << "CollectiveFloatFun '" << getName() << "': Error! cutiter, floatiter and intiter all 0!\n";
    }
  }
  
  if (!intiter) {
    std::cerr << "CollectiveFloatFun '" << getName() << "': Cannot determine intiter from fun and cut!\n";
  }
  
  assert (intiter);
  assert (!extiter || extiter!=intiter);
  assert (!funiter || funiter==intiter || funiter==extiter);
  assert (!cutiter || cutiter==intiter || cutiter==extiter);
}

CachedFloatFun::CachedFloatFun (const ROListPoR& rol,  
                                const FloatFun&  arg_ 
                               )
: FloatFun (std::string("cached(")+arg_.getName()+')'),
  CachedO (rol), 
  arg(arg_), 
  cacheValid (false)                 
{ 
   assert (&arg); 
}   

Float_FF CachedFloatFun::operator() () const {
  if (!cacheValid) {
    cachedValue = arg();
    cacheValid = true;
  }
  return cachedValue;
}

CachedFloatFunIt::CachedFloatFunIt (const ROListPoR&     rol, 
                                    const FloatFun&      arg_, 
                                    const FillIterator&  iter_ 
                                   )
: FloatFun (std::string("cached(")+arg_.getName()+','+iter_.getName()+')'),
  CachedO (rol), 
  arg (arg_), 
  iter (iter_), 
  cacheValid (), 
  cachedValues ()              
{ 
  assert (&arg); 
  assert (&iter);
}   

Float_FF CachedFloatFunIt::operator() () const {
  int i = iter();
  assert (i >= 0);
  if (static_cast<unsigned int>(i) >= cacheValid.size()) growCache (i+1); 
  if (!cacheValid[i]) {
    cacheValid[i] = true;
    cachedValues[i] = arg();
  }
  return cachedValues[i];
}

void CachedFloatFunIt::growCache (unsigned int s) const {
  unsigned int oldsize = cacheValid.size();
  if (s < oldsize) return;
  cacheValid.resize (s);
  cachedValues.resize(s);
  for (unsigned int i = oldsize; i < s; ++i) cacheValid[i] = false;
}

FloatFun& cached (const ROListPoR&   rol,    
                  const FloatFunPoR& fun,    
                  const FillIteratorPoR& iter
                 ) {
  assert (fun.pff);                     
  if (iter.pfi) return *new CachedFloatFunIt (rol, *fun.pff, *iter.pfi);
  if (const FillIterator *fi = fun.pff->getIterator()) 
    return  *new CachedFloatFunIt (rol, *fun.pff, *fi);
  return *new CachedFloatFun (rol, *fun.pff);                  
}

FloatFun& operator+ (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new SumOfTwoFuns (lhs_, rhs_);
}

FloatFun& operator+ (const FloatFun& lhs_, Float_FF rhs_) {
  return *new ScalarSumFloatFun (rhs_, lhs_);
}

FloatFun& operator+ (Float_FF lhs_, const FloatFun& rhs_) {
  return *new ScalarSumFloatFun (lhs_, rhs_);
}

FloatFun& operator- (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new DiffOfTwoFuns (lhs_, rhs_);
}

FloatFun& operator- (const FloatFun& lhs_, Float_FF rhs_) {
  return *new ScalarSumFloatFun (-rhs_, lhs_, 1, NamedO::str('(')+lhs_.getName()+")-"+NamedO::str(rhs_));
}

FloatFun& operator- (Float_FF lhs_, const FloatFun& rhs_) {
  return *new ScalarSumFloatFun (lhs_, rhs_, -1, NamedO::str(lhs_)+"-("+rhs_.getName()+')');
}

FloatFun& operator- (const FloatFun& arg_) {
  return *new ScalarProductFloatFun ("-", -1, arg_);
}

const FloatFun& operator+ (const FloatFun& arg_) {
  return arg_;
}
FloatFun& operator+ (FloatFun& arg_) {
  return arg_;
}

FloatFun& operator* (const FloatFun& lhs_, const FloatFun& rhs_) {
  if (&lhs_ == &rhs_) return square (lhs_);
  return *new ProdOfTwoFuns (lhs_, rhs_);
}

FloatFun& operator* (Float_FF lhs_, const FloatFun& rhs_) {
  return *new ScalarProductFloatFun (lhs_, rhs_);
}

FloatFun& operator* (const FloatFun& lhs_, Float_FF rhs_) {
  return *new ScalarProductFloatFun (rhs_, lhs_);
}

FloatFun& operator/ (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new DivOfTwoFuns (lhs_, rhs_);
}

FloatFun& operator/ (Float_FF lhs_, const FloatFun& rhs_) {
  return *new ScalarDivisionFloatFun (lhs_, rhs_);
}

FloatFun& operator/ (const FloatFun& lhs_, Float_FF rhs_) {
  return *new ScalarProductFloatFun (1./rhs_, lhs_);
}



FloatFunComp& operator< (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunLess (lhs_, rhs_);
}

FloatFunComp& operator> (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunGreater (lhs_, rhs_);
}

FloatFunComp& operator<= (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunLessEqual (lhs_, rhs_);
}

FloatFunComp& operator>= (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunGreaterEqual (lhs_, rhs_);
}

FloatFunComp& operator== (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunEqual (lhs_, rhs_);
}

FloatFunComp& operator!= (const FloatFun& lhs_, const FloatFun& rhs_) {
  return *new FloatFunNotEqual (rhs_, lhs_);
}

// External functions that implement comparison operators between
// a Float_FF and a FloatFun object

FloatFunComp& operator< (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunLess (lhs_, rhs_);
}

FloatFunComp& operator<= (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunLessEqual (lhs_, rhs_);
}

FloatFunComp& operator> (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunGreater (lhs_, rhs_);
}

FloatFunComp& operator>= (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunGreaterEqual (lhs_, rhs_);
}

FloatFunComp& operator== (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunEqual (lhs_, rhs_);
}

FloatFunComp& operator!= (Float_FF lhs_, const FloatFun& rhs_) {
  return *new FloatFloatFunNotEqual (lhs_, rhs_);
}

// External functions that implement comparison operators between
// a FloatFun object and a Float_FF

FloatFunComp& operator< (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatLess (lhs_, rhs_);
}

FloatFunComp& operator<= (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatLessEqual (lhs_, rhs_);
}

FloatFunComp& operator> (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatGreater (lhs_, rhs_);
}

FloatFunComp& operator>= (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatGreaterEqual (lhs_, rhs_);
}

FloatFunComp& operator== (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatEqual (lhs_, rhs_);
}

FloatFunComp& operator!= (const FloatFun& lhs_, Float_FF rhs_) {
  return *new FloatFunFloatNotEqual (lhs_, rhs_);
}

// External functions that allow chain comparisons like 1 < x < y <= 15
// Such a chain comparison is defined as
// (1 < x) && (x < y) && (y <= 15)

FloatFunComp& operator< (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompLess (lhs_, rhs_);
}

FloatFunComp& operator<= (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompLessEqual (lhs_, rhs_);
}

FloatFunComp& operator> (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompGreater (lhs_, rhs_);
}

FloatFunComp& operator>= (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompGreaterEqual (lhs_, rhs_);
}

FloatFunComp& operator== (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompEqual (lhs_, rhs_);
}

FloatFunComp& operator!= (const FloatFunComp& lhs_, const FloatFun& rhs_) {
  return *new ChainCompNotEqual (lhs_, rhs_);
}

FloatFunComp& operator< (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatLess (lhs_, rhs_);
}

FloatFunComp& operator<= (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatLessEqual (lhs_, rhs_);
}

FloatFunComp& operator> (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatGreater (lhs_, rhs_);
}

FloatFunComp& operator>= (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatGreaterEqual (lhs_, rhs_);
}

FloatFunComp& operator== (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatEqual (lhs_, rhs_);
}

FloatFunComp& operator!= (const FloatFunComp& lhs_, Float_FF rhs_) {
  return *new ChainCompFloatNotEqual (lhs_, rhs_);
}


FloatFun& max (const FloatFunPoRConst&      fun,      
               const BaseCutPoRConst&       cut,      
               const FillIteratorPoR&       intiter,  
               const FillIteratorPoRConst&  extiter   
              ) {
  return *new MaxFloatFun (fun, cut, intiter, extiter);
}

FloatFun& min (const FloatFunPoRConst&      fun,        
               const BaseCutPoRConst&       cut,     
               const FillIteratorPoR&       intiter, 
               const FillIteratorPoRConst&  extiter  
              ) {
  return *new MinFloatFun (fun, cut, intiter, extiter);
}

FloatFun& sum (const FloatFunPoRConst&      fun,      
               const BaseCutPoRConst&       cut,     
               const FillIteratorPoR&       intiter, 
               const FillIteratorPoRConst&  extiter  
              ) {
  return *new SumFloatFun (fun, cut, intiter, extiter);
}

FloatFun& average (const FloatFunPoRConst&      fun,       
                   const BaseCutPoRConst&       cut,     
                   const FillIteratorPoR&       intiter, 
                   const FillIteratorPoRConst&  extiter  
                  ) {
  return *new AverageFloatFun (fun, cut, intiter, extiter);
}

---