CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/JetMETCorrections/InterpolationTables/interface/StorableInterpolationFunctor.h

Go to the documentation of this file.
00001 #ifndef NPSTAT_STORABLEINTERPOLATIONFUNCTOR_HH_
00002 #define NPSTAT_STORABLEINTERPOLATIONFUNCTOR_HH_
00003 
00015 #include "JetMETCorrections/InterpolationTables/interface/StorableMultivariateFunctor.h"
00016 #include "JetMETCorrections/InterpolationTables/interface/LinInterpolatedTableND.h"
00017 #include "JetMETCorrections/InterpolationTables/interface/SimpleFunctors.h"
00018 
00019 namespace npstat {
00024     template
00025     <
00026         class Numeric,
00027         class Axis = UniformAxis,
00028         class Converter = Same<Numeric>
00029     >
00030     class StorableInterpolationFunctor : public StorableMultivariateFunctor
00031     {
00032         template <typename Num2, typename Axis2, typename Conv2>
00033         friend class StorableInterpolationFunctor;
00034 
00035     public:
00036         typedef LinInterpolatedTableND<Numeric,Axis> Table;
00037 
00039 
00040         template <class Num2>
00041         inline StorableInterpolationFunctor(
00042             const LinInterpolatedTableND<Num2,Axis>& table)
00043             : StorableMultivariateFunctor(), table_(table) {}
00044 
00045         template <class Num2>
00046         inline StorableInterpolationFunctor(
00047             const LinInterpolatedTableND<Num2,Axis>& table,
00048             const std::string& descr)
00049             : StorableMultivariateFunctor(descr), table_(table) {}
00051 
00053         template <class Num2, class Conv2>
00054         inline StorableInterpolationFunctor(
00055             const StorableInterpolationFunctor<Num2,Axis,Conv2>& tab)
00056             : StorableMultivariateFunctor(tab.description()),
00057               table_(tab.table_) {}
00058 
00060 
00065         inline StorableInterpolationFunctor(
00066             const std::vector<Axis>& axes,
00067             const std::vector<std::pair<bool,bool> >& interpolationType,
00068             const char* functionLabel=0)
00069             : StorableMultivariateFunctor(),
00070               table_(axes, interpolationType, functionLabel) {}
00071 
00072         inline StorableInterpolationFunctor(
00073             const Axis& xAxis, bool leftX, bool rightX,
00074             const char* functionLabel=0)
00075             : StorableMultivariateFunctor(),
00076               table_(xAxis, leftX, rightX, functionLabel) {}
00077 
00078         inline StorableInterpolationFunctor(
00079             const Axis& xAxis, bool leftX, bool rightX,
00080             const Axis& yAxis, bool leftY, bool rightY,
00081             const char* functionLabel=0)
00082             : StorableMultivariateFunctor(),
00083               table_(xAxis, leftX, rightX,
00084                      yAxis, leftY, rightY, functionLabel) {}
00085 
00086         inline StorableInterpolationFunctor(
00087             const Axis& xAxis, bool leftX, bool rightX,
00088             const Axis& yAxis, bool leftY, bool rightY,
00089             const Axis& zAxis, bool leftZ, bool rightZ,
00090             const char* functionLabel=0)
00091             : StorableMultivariateFunctor(),
00092               table_(xAxis, leftX, rightX,
00093                      yAxis, leftY, rightY,
00094                      zAxis, leftZ, rightZ, functionLabel) {}
00095 
00096         inline StorableInterpolationFunctor(
00097             const Axis& xAxis, bool leftX, bool rightX,
00098             const Axis& yAxis, bool leftY, bool rightY,
00099             const Axis& zAxis, bool leftZ, bool rightZ,
00100             const Axis& tAxis, bool leftT, bool rightT,
00101             const char* functionLabel=0)
00102             : StorableMultivariateFunctor(),
00103               table_(xAxis, leftX, rightX,
00104                      yAxis, leftY, rightY,
00105                      zAxis, leftZ, rightZ,
00106                      tAxis, leftT, rightT, functionLabel) {}
00107 
00108         inline StorableInterpolationFunctor(
00109             const Axis& xAxis, bool leftX, bool rightX,
00110             const Axis& yAxis, bool leftY, bool rightY,
00111             const Axis& zAxis, bool leftZ, bool rightZ,
00112             const Axis& tAxis, bool leftT, bool rightT,
00113             const Axis& vAxis, bool leftV, bool rightV,
00114             const char* functionLabel=0)
00115             : StorableMultivariateFunctor(),
00116               table_(xAxis, leftX, rightX,
00117                      yAxis, leftY, rightY,
00118                      zAxis, leftZ, rightZ,
00119                      tAxis, leftT, rightT,
00120                      vAxis, leftV, rightV, functionLabel) {}
00122 
00123         virtual ~StorableInterpolationFunctor() {}
00124 
00125         virtual unsigned minDim() const {return table_.dim();};
00126 
00127         virtual double operator()(const double* point, unsigned dim) const
00128             {return conv_(table_(point, dim));}
00129 
00131 
00132         inline Table& interpolator() {return table_;}
00133         inline const Table& interpolator() const {return table_;}
00135 
00137 
00138         inline ArrayND<Numeric>& table() {return table_.table();}
00139         inline const ArrayND<Numeric>& table() const {return table_.table();}
00141 
00143         inline void setConverter(const Converter& conv) {conv_ = conv;}
00144 
00146         // Method related to "geners" I/O
00147         virtual gs::ClassId classId() const {return gs::ClassId(*this);}
00148         virtual bool write(std::ostream& of) const;
00150 
00151         // I/O methods needed for reading
00152         static inline const char* classname();
00153         static inline unsigned version() {return 1;}
00154         static StorableInterpolationFunctor* read(
00155             const gs::ClassId& id, std::istream& in);
00156 
00157     protected:
00158         virtual bool isEqual(const StorableMultivariateFunctor& other) const
00159         {
00160             // Note the use of static_cast rather than dynamic_cast below.
00161             // static_cast works faster and it is guaranteed to succeed here.
00162             const StorableInterpolationFunctor& r = 
00163                 static_cast<const StorableInterpolationFunctor&>(other);
00164             return table_ == r.table_ &&
00165                    this->description() == other.description();
00166         }
00167 
00168     private:
00169         StorableInterpolationFunctor();
00170 
00171         Table table_;
00172         Converter conv_;
00173     };
00174 }
00175 
00176 #include "Alignment/Geners/interface/binaryIO.hh"
00177 #include "Alignment/Geners/interface/CPP11_auto_ptr.hh"
00178 #include "Alignment/Geners/interface/IOException.hh"
00179 
00180 namespace npstat {
00181     template <typename Numeric, class Axis, class Converter>
00182     const char* StorableInterpolationFunctor<Numeric,Axis,Converter>::classname()
00183     {
00184         static const std::string myClass(gs::template_class_name<Numeric,Axis>(
00185                                       "npstat::StorableInterpolationFunctor"));
00186         return myClass.c_str();
00187     }    
00188 
00189     template<typename Numeric, class Axis, class Converter>
00190     bool StorableInterpolationFunctor<Numeric,Axis,Converter>::write(
00191         std::ostream& of) const
00192     {
00193         gs::write_pod(of, this->description());
00194         return table_.classId().write(of) && table_.write(of);
00195     }
00196 
00197     template<typename Numeric, class Axis, class Converter>
00198     StorableInterpolationFunctor<Numeric,Axis,Converter>*
00199     StorableInterpolationFunctor<Numeric,Axis,Converter>::read(
00200         const gs::ClassId& id, std::istream& in)
00201     {
00202         static const gs::ClassId current(
00203             gs::ClassId::makeId<StorableInterpolationFunctor<Numeric,Axis> >());
00204         current.ensureSameId(id);
00205 
00206         std::string descr;
00207         gs::read_pod(in, &descr);
00208         gs::ClassId tabid(in, 1);
00209         if (in.fail()) throw gs::IOReadFailure(
00210             "In npstat::StorableInterpolationFunctor::read: "
00211             "input stream failure");
00212         CPP11_auto_ptr<Table> tab(Table::read(tabid, in));
00213         return new StorableInterpolationFunctor(*tab, descr);
00214     }
00215 }
00216 
00217 
00218 #endif // NPSTAT_STORABLEINTERPOLATIONFUNCTOR_HH_
00219