ROOT logo
#ifndef  __HPARTICLECUT_H__
#define  __HPARTICLECUT_H__

#include "TString.h"
#include "TNamed.h"
#include "TTree.h"
#include "TTreeFormula.h"
#include "TDirectory.h"
#include "TROOT.h"


#include <iostream>
#include <iomanip>
#include <vector>

using namespace std;


template <class T>
class HParticleCut : public TNamed {

private:
    Int_t         fCutNumber;    // a unique cut number to handle this cut
    TString       fcondition;    // formular expression for this object to cut on
    UInt_t           fmaxCut;    // how many version can be used at max
    TDirectory*     fsaveDir;    //!  remember the actual directory before strting to work (will be restored)
    vector<ULong64_t>fctFail;    // count evaluation == kFALSE
    vector<ULong64_t>fctCall;    // count all call to the evaluation
    Bool_t      fbInverseCut;    // default kFALSE , kTRUE will invert selection
    T*                    fc;    // template class pointer
    TTree*             fTree;    // mini tree on the object to cut on
    TTreeFormula*    fselect;    // formula object

public:


    HParticleCut(TString name,Int_t num,TString cond);
    ~HParticleCut();
    void          setInverse(Bool_t inv)   { fbInverseCut = inv;}
    void          setMaxCut(UInt_t max=4);
    Bool_t        eval(T* c,UInt_t version=0);
    void          print();
    void          resetCounter();
    Float_t       getCutRate(UInt_t version=0);
    UInt_t        getNCall(UInt_t version=0);
    UInt_t        getNFail(UInt_t version=0);
    T*            getClass()     { return fc; }
    TTree*        getTree()      { return fTree; }
    TTreeFormula* getFormula()   { return fselect;}
    TString       getCondition() { return fcondition;}
    ClassDef(HParticleCut,1)
};

template <class T>
HParticleCut<T>::HParticleCut(TString name,Int_t num,TString cond)
{
    fsaveDir = gDirectory;
    gROOT->cd();
    SetName(name);
    fCutNumber = num;
    fcondition = cond;
    fbInverseCut = kFALSE;
    fc = new T();
    fTree = new TTree(Form("Tree_%s_%i",name.Data(),fCutNumber),Form("Tree_%s_%i",name.Data(),fCutNumber));
    fTree->Branch(fc->ClassName(),fc->ClassName(),&fc,99);
    fselect = new TTreeFormula(Form("Selection_%s_%i",name.Data(),fCutNumber),fcondition,fTree);
    setMaxCut(4);
    fsaveDir->cd();

}

template <class T>
HParticleCut<T>::~HParticleCut()
{
    if(fselect)delete fselect;
    if(fTree)  delete fTree;
    if(fc)     delete fc;
}


template <class T>
void   HParticleCut<T>::setMaxCut(UInt_t max)
{
    // set the max number of versions of stst counters
    fmaxCut     = max;
    if(fmaxCut == 0) fmaxCut = 1;
    resetCounter();
}


template <class T>
Float_t HParticleCut<T>::getCutRate(UInt_t version)
{
    return (version <fmaxCut && fctCall[version] > 0) ? (fctFail[version]/(Float_t)fctCall[version])*100.:0;
}

template <class T>
UInt_t HParticleCut<T>::getNCall(UInt_t version)
{
    return (version <fmaxCut ) ? fctCall[version]:0;
}

template <class T>
UInt_t HParticleCut<T>::getNFail(UInt_t version)
{
    return (version <fmaxCut ) ? fctFail[version]:0;
}

template <class T>
void HParticleCut<T>::resetCounter()
{
    fctFail.clear();
    fctCall.clear();
    fctFail.assign(fmaxCut,0);
    fctCall.assign(fmaxCut,0);
}

template <class T>
Bool_t HParticleCut<T>::eval(T* c,UInt_t version)
{
    // returns kTRUE if the object fullfills the condition
    // counters for calls and failed conditions are filled
    if(version<fmaxCut){
	fctCall[version]++;
    }

    new (fc) T(*c);

    if ( (!fbInverseCut && fselect->EvalInstance(0) == 0) ||
         ( fbInverseCut && fselect->EvalInstance(0) != 0)
       )
    {
	if(version<fmaxCut)fctFail[version]++;
	return kFALSE;

    } else   { return kTRUE; }
}

template <class T>
void HParticleCut<T>::print()
{
    Int_t p = cout.precision();
    std::ios_base::fmtflags fl =cout.flags();

    cout<<"CutNumber : "<<setw(3)<<fCutNumber
	<<" name : "<<setw(20)<<GetName()
	<<", cut [%] : "<<flush;
    for(UInt_t i = 0; i < fmaxCut ; i++){ cout<<setw(5)<<fixed<<right<<setprecision(1)<<getCutRate(i)<<" "<<flush;}
    cout<<", cut = "<<fcondition <<setprecision(p)<<endl;
    cout.flags(fl);

}

#endif
 hparticlecut.h:1
 hparticlecut.h:2
 hparticlecut.h:3
 hparticlecut.h:4
 hparticlecut.h:5
 hparticlecut.h:6
 hparticlecut.h:7
 hparticlecut.h:8
 hparticlecut.h:9
 hparticlecut.h:10
 hparticlecut.h:11
 hparticlecut.h:12
 hparticlecut.h:13
 hparticlecut.h:14
 hparticlecut.h:15
 hparticlecut.h:16
 hparticlecut.h:17
 hparticlecut.h:18
 hparticlecut.h:19
 hparticlecut.h:20
 hparticlecut.h:21
 hparticlecut.h:22
 hparticlecut.h:23
 hparticlecut.h:24
 hparticlecut.h:25
 hparticlecut.h:26
 hparticlecut.h:27
 hparticlecut.h:28
 hparticlecut.h:29
 hparticlecut.h:30
 hparticlecut.h:31
 hparticlecut.h:32
 hparticlecut.h:33
 hparticlecut.h:34
 hparticlecut.h:35
 hparticlecut.h:36
 hparticlecut.h:37
 hparticlecut.h:38
 hparticlecut.h:39
 hparticlecut.h:40
 hparticlecut.h:41
 hparticlecut.h:42
 hparticlecut.h:43
 hparticlecut.h:44
 hparticlecut.h:45
 hparticlecut.h:46
 hparticlecut.h:47
 hparticlecut.h:48
 hparticlecut.h:49
 hparticlecut.h:50
 hparticlecut.h:51
 hparticlecut.h:52
 hparticlecut.h:53
 hparticlecut.h:54
 hparticlecut.h:55
 hparticlecut.h:56
 hparticlecut.h:57
 hparticlecut.h:58
 hparticlecut.h:59
 hparticlecut.h:60
 hparticlecut.h:61
 hparticlecut.h:62
 hparticlecut.h:63
 hparticlecut.h:64
 hparticlecut.h:65
 hparticlecut.h:66
 hparticlecut.h:67
 hparticlecut.h:68
 hparticlecut.h:69
 hparticlecut.h:70
 hparticlecut.h:71
 hparticlecut.h:72
 hparticlecut.h:73
 hparticlecut.h:74
 hparticlecut.h:75
 hparticlecut.h:76
 hparticlecut.h:77
 hparticlecut.h:78
 hparticlecut.h:79
 hparticlecut.h:80
 hparticlecut.h:81
 hparticlecut.h:82
 hparticlecut.h:83
 hparticlecut.h:84
 hparticlecut.h:85
 hparticlecut.h:86
 hparticlecut.h:87
 hparticlecut.h:88
 hparticlecut.h:89
 hparticlecut.h:90
 hparticlecut.h:91
 hparticlecut.h:92
 hparticlecut.h:93
 hparticlecut.h:94
 hparticlecut.h:95
 hparticlecut.h:96
 hparticlecut.h:97
 hparticlecut.h:98
 hparticlecut.h:99
 hparticlecut.h:100
 hparticlecut.h:101
 hparticlecut.h:102
 hparticlecut.h:103
 hparticlecut.h:104
 hparticlecut.h:105
 hparticlecut.h:106
 hparticlecut.h:107
 hparticlecut.h:108
 hparticlecut.h:109
 hparticlecut.h:110
 hparticlecut.h:111
 hparticlecut.h:112
 hparticlecut.h:113
 hparticlecut.h:114
 hparticlecut.h:115
 hparticlecut.h:116
 hparticlecut.h:117
 hparticlecut.h:118
 hparticlecut.h:119
 hparticlecut.h:120
 hparticlecut.h:121
 hparticlecut.h:122
 hparticlecut.h:123
 hparticlecut.h:124
 hparticlecut.h:125
 hparticlecut.h:126
 hparticlecut.h:127
 hparticlecut.h:128
 hparticlecut.h:129
 hparticlecut.h:130
 hparticlecut.h:131
 hparticlecut.h:132
 hparticlecut.h:133
 hparticlecut.h:134
 hparticlecut.h:135
 hparticlecut.h:136
 hparticlecut.h:137
 hparticlecut.h:138
 hparticlecut.h:139
 hparticlecut.h:140
 hparticlecut.h:141
 hparticlecut.h:142
 hparticlecut.h:143
 hparticlecut.h:144
 hparticlecut.h:145
 hparticlecut.h:146
 hparticlecut.h:147
 hparticlecut.h:148
 hparticlecut.h:149
 hparticlecut.h:150
 hparticlecut.h:151
 hparticlecut.h:152
 hparticlecut.h:153
 hparticlecut.h:154