ROOT logo
HYDRA - THE HADES ANALYSIS PACKAGE » (UNKNOWN) » HParticleEvtChara

class HParticleEvtChara: public HReconstructor

_HADES_CLASS_DESCRIPTION 
*-- AUTHOR : B. Kardan  28.08.2019
*-- VERSION : 0.73

 HParticleEvtChara

 Purpose: EventCharacterization
 - Centrality from Hit(TOF,RPC,FW) and Track Estimators in Data/Sim
 - QVector and Phi (Reaction Plane Estimate) from FW
 - Event-weight for downscaled Events PT2


 Usage:

  - input files can be found at : example /cvmfs/hades.gsi.de/param/eventchara

  - to define the ParameterFile where Classes and Estimators are stored use
    setParameterFile("/cvmfs/hades.gsi.de/param/eventchara/ParameterFile.root")

  - to print the definition of estimator & class use
    printCentralityClass(HParticleEvtChara::kTOFRPC, HParticleEvtChara::k10);

  - to get the CentralityClass of an event (with estimator and class definition) use
    getCentralityClass(HParticleEvtChara::kTOFRPC, HParticleEvtChara::k10);

  - to get the CentralityPercentile of an event (only estimator is needed) use
    getCentralityPercentile(HParticleEvtChara::kTOFRPC);

  - to get the EventWeight to re-weight downscaled events use
    getEventWeight();

  - to get the EventPlane with ReCentering use
    getEventPlane(HParticleEvtChara::kDefault);


   Estimators:
    TOFRPC                   - (default) TOF and RPC hit multiplicity in time-window
    TOF                      - TOF hit multiplicity in time-window
    RPC                      - RPC hit multiplicity in time-window
    TOFRPCtot                - total TOF and RPC hit multiplicity in event-window
    SelectedParticleCand     - selected Particle multiplicity
    PrimaryParticleCand      - primary Particle multiplicity
    SelectedParticleCandCorr - selected Particle multiplicity corrected by the
                               running mean and scaled with the reference mean
                               (selTrack * referenceMean/<selTrack>)
    SelectedParticleCandSecCorr
                              - selected Particle multiplicity corrected by the
                               running mean and scaled with the reference mean
                               (selTrack * referenceMean/<selTrack>)
    SelectedParticleCandNorm - selected Particle multiplicity normalized by
                               the running mean (selTrack/<selTrack>)
    FWSumChargeSpec          - sum of charge (dE/dx in a.u.) of Spectator in
                               FW acceptance with beta~0.9
    FWSumChargeZ             - sum of charge (int Z up till charge-state 14
                               with individuel fixed cuts in dE/dx each FW-cell)
                               of Spectator in FW acceptance with beta~0.9

   Classes:
    2                   - 2% classes
    5                   - 5% classes
    10                  - (default) 10% classes
    13                  - 13% classes
    15                  - 15% classes
    20                  - 20% classes
    FOPI                - FOPI centrality classes

  EventPlaneCorrection:
    kNoCorrection            - EP only selection of spectator candidates in FW,
                               no Correction
    kSmear                   - smearing of x and y of each FW hit inside cell size
    kShiftFW                 - global shift x and y with the centre of gravity and
                               scale x and y with the sigma of the distribution of
                               all FW hits in all events (per class/day)
    kWeightCharge            - weigthing with charge-state up to 14 with individuel
                               fixed timing- and dE/dx-cuts each FW-cell
    kReCentering             - re-centering of QVector with <Qx> and <Qy> (calc. evt-by-evt)
                               (only first harmonic correction)
    kScaling                 - scaling of QVector with the sigma of Qx and Qy (calc. evt-by-evt)
    kRotation                - rotation of EP via residual Fourier harmonics up to 8 cos and sin terms
                               after FWshift and FWscaling

    kDefault                 - recommended option (kShiftFW|kWeightCharge|kRotation)

 quick how-To:
  - for an example see demo-macro:
       /scripts/batch/GE/loopDST/loopDST.C

  1. include header-file
      #include "HParticleEvtChara.h"

  2. if you use Hloop check that following catagories are loaded:

   HParticleEvtInfo
   HParticleCand
   HWallHit

       if(!loop.setInput("-*,+HParticleEvtInfo, +HParticleCand, +HWallHit")){
           cout<<"READBACK: ERROR : cannot read input !"<<endl;
       }

  3. before eventLoop:

       HParticleEvtChara evtChara;
       TString ParameterfileCVMFS = "/lustre/nyx/hades/user/bkardan/param/centrality_epcorr_apr12_gen8_2019_02_pass30.root";

       //due to the overlap of day126  there is dedidacate param-file for the reverse-field runs
       if(isRevFieldDATA)   ParameterfileCVMFS = "/lustre/nyx/hades/user/bkardan/param/centrality_epcorr_apr12_gen8_revfield_2019_02_pass29.root";

       //if you need centrality-values for sim UrQMD gen8a you find them here
       //due missing primary light nuclei in UrQMD... there is no correction on the FW-EventPlane in this version!
       if(isSimulation)   ParameterfileCVMFS = "/lustre/nyx/hades/user/bkardan/param/centrality_epcorr_sim_au1230au_gen8a_UrQMD_minbias_2019_04_pass0.root";
       //if(isSimulation) ParameterfileCVMFS = "/lustre/nyx/hades/user/bkardan/param/centrality_epcorr_sim_au1230au_gen9vertex_UrQMD_minbias_2019_04_pass0.root";

       if(!evtChara.setParameterFile(ParameterfileCVMFS)){
           cout << "Parameterfile not found !!! " << endl;return kFALSE;
       }

       if(!evtChara.init()) {
           cout << "HParticleEvtChara not init!!! " << endl;return kFALSE;
       }

       Int_t eCentEst   = HParticleEvtChara::kTOFRPC;
       Int_t eCentClass = HParticleEvtChara::k10;
       Int_t eEPcorr    = HParticleEvtChara::kDefault;
       cout << "\t selected EPcorrection method is:  "  << evtChara.getStringEventPlaneCorrection(eEPcorr) << endl;


       evtChara.printCentralityClass(eCentEst, eCentClass);

  4. inside event-loop:

        Float_t  event_weight = evtChara.getEventWeight();   // event_weight dependent if PT2(down-scaled) or PT3

        Int_t   CentralityClassTOFRPC = evtChara.getCentralityClass(eCentEst, eCentClass);
        // 10% Centrality-Classes:       1(0-10%) - 5(40-50%) ... 0:Overflow max:Underflow

        Float_t CentralityTOFRPC      = evtChara.getCentralityPercentile(eCentClass);
        // CentralityPercentile:         0 - 100% percentile of the total cross section
        //                                   101% Over-,Underflow or Outlier Events

        Float_t EventPlane            = evtChara.getEventPlane(eEPcorr);
        //EventPlane:                   0 - 2pi (in rad)  re-centered & scaled EP
        //                             -1   no EP could be reconstructed

        //for EP-resolution use the sub-event:
        Float_t EventPlaneA            = evtChara.getEventPlane(eEPcorr,1);
        Float_t EventPlaneB            = evtChara.getEventPlane(eEPcorr,2);

        //check if the EP(and subEvent EP) could be reconstructed, most of the case less than 4 Hits in FW
        if(EventPlane  == -1)  continue;
        if(EventPlaneA  == -1 || EventPlaneB  == -1 ) continue;

        check if you have any further selections on FW-hit multiplicites in you own code,
        this will bias the result!

  5. to calculate the EP-resolution one option is the following:
     (Ollitrault method.  - approximation valid for high chi-values)
     inside event-loop fill the difference of the two sub-event EP
     into a histogramm for each centrality class:

          Float_t deltaPhiAB =  TMath::Abs(TVector2::Phi_mpi_pi(EventPlaneA-EventPlaneB));
          h->Fill(deltaPhiAB);

     after your event-analysis you can finalize your results by calculating the EP-resolution for each centrality class
     by the ratio of event N(pi/2 < deltaEP_AB < pi)/N(0 < deltaEP_AB < pi):

         Int_t i90        = h->FindBin(0.5*TMath::Pi());
         Int_t i180       = h->FindBin(TMath::Pi());
         Double_t ratio   = h->Integral(i90,i180) /h->Integral(1,i180);
         Double_t dChi    = sqrt(-2.*TMath::Log(2.*ratio));
         Double_t dEPres1 = HParticleEvtChara::ComputeResolution( dChi );    // for first harmonic
         Double_t dEPres2 = HParticleEvtChara::ComputeResolution( dChi , 2); // for second ...

  6. the second option will be described here soon here ;)
     (Voloshin method - precise iterative resolution search)
     inside event-loop fill the Cosine of the difference of the two sub-event EP
     into a histogramm for each centrality class:

         Float_t CosDeltaPhiAB =  TMath::Cos(EventPlaneA-EventPlaneB);
         h2->Fill(CosDeltaPhiAB);

     after your event-analysis you can finalize your results by calculating the EP-resolution for each centrality class
     by the mean-value:

         Double_t MeanCosDeltaPhiAB  = h2->GetMean();        //< cos n delta_phi>;
         Double_t dV = TMath::Sqrt(MeanCosDeltaPhiAB);
         Double_t dChi = FindXi(dV,1e-6, 1);                 // for first harmonic
         Double_t dEPres1 = HParticleEvtChara::ComputeResolution( dChi );    // for first harmonic
         Double_t dEPres2 = HParticleEvtChara::ComputeResolution( dChi , 2); // for second ...
         dV = ComputeResolution( TMath::Sqrt2()*dChi , 1);
         printf("An estimate of the event plane resolution first order is: %f\n", dV );


  7. to correct your flow-values use 1./dEPres1

  8. test your macro with data and if you real have the right parameterfile
      with right values for centrality and EPcorections, check the output of

        evtChara.printCentralityClass(eCentEst, eCentClass);

  9. and now happy plotting ;)


 Change History:

 19.02.2019  release of 0.7
 24.04.2019  implementation of useFWCut() and check if FWcut hist are loaded
 04.06.2019  fix of embeded-sim into data and removed of high-momentum condition

Function Members (Methods)

public:
HParticleEvtChara(const Text_t* name = "HParticleEvtChara", const Text_t* title = "HParticleEvtChara")
virtual~HParticleEvtChara()
voidTObject::AbstractMethod(const char* method) const
Bool_taddEstimatorHist(TH1F* hist, Float_t fractionXsection = 100., UInt_t centE = kTOFRPC, Int_t direction = -1)
Bool_taddEventPlaneCorrectionHist(TProfile2D* hist, UInt_t epParam = kQx)
Bool_taddFWCutValuesHist(TH1* hist, Int_t cell, UInt_t eFWCut)
virtual voidTObject::AppendPad(Option_t* option = "")
virtual voidHReconstructor::Browse(TBrowser* b)
static TClass*Class()
virtual const char*TObject::ClassName() const
virtual voidTNamed::Clear(Option_t* option = "")
virtual TObject*TNamed::Clone(const char* newname = "") const
virtual Int_tTNamed::Compare(const TObject* obj) const
Double_tComputeResolution(Double_t x, Int_t n = 1) const
virtual Bool_tHReconstructor::connectTask(HTask* task, Int_t n = 0)
virtual voidTNamed::Copy(TObject& named) const
virtual voidTObject::Delete(Option_t* option = "")MENU
virtual Int_tTObject::DistancetoPrimitive(Int_t px, Int_t py)
virtual voidTObject::Draw(Option_t* option = "")
virtual voidTObject::DrawClass() constMENU
virtual TObject*TObject::DrawClone(Option_t* option = "") constMENU
virtual voidTObject::Dump() constMENU
virtual voidTObject::Error(const char* method, const char* msgfmt) const
virtual Int_texecute()
virtual voidTObject::Execute(const char* method, const char* params, Int_t* error = 0)
virtual voidTObject::Execute(TMethod* method, TObjArray* params, Int_t* error = 0)
virtual voidTObject::ExecuteEvent(Int_t event, Int_t px, Int_t py)
virtual voidTObject::Fatal(const char* method, const char* msgfmt) const
virtual voidTNamed::FillBuffer(char*& buffer)
voidfillHitArray()
Bool_tfillQVectors(UInt_t EPcorr = kNoCorrection, UInt_t nHarmonic = 1)
virtual Bool_tfinalize()
virtual TObject*TObject::FindObject(const char* name) const
virtual TObject*TObject::FindObject(const TObject* obj) const
Double_tFindXi(Double_t res, Double_t prec = 1e-6, Int_t n = 1) const
Float_t*getBinCenterArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Int_tgetCentralityClass(TString estimator)
Int_tgetCentralityClass(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Int_tgetCentralityClass10(TString estimator)
Int_tgetCentralityClass5(TString estimator)
Float_t*getCentralityClassArray(UInt_t centC = k10)
Float_tgetCentralityClassBinSize(UInt_t e)
TH1F*getCentralityClassHist(UInt_t centE = kTOFRPC, UInt_t centC = k10) const
Int_tgetCentralityClassNbins(UInt_t centC = k10)
Float_tgetCentralityEstimator(UInt_t centE = kTOFRPC)
Int_tgetCentralityEstimatorBinMax(UInt_t e)
Int_tgetCentralityEstimatorBinSize(UInt_t e)
UInt_tgetCentralityEstimatorFromString(TString s)
Float_tgetCentralityPercentile(TString estimator)
Float_tgetCentralityPercentile(UInt_t centE = kTOFRPC)
TH1F*getCentralityPercentileHist(UInt_t centE = kTOFRPC) const
virtual HTask*HReconstructor::getComposite()
virtual voidHReconstructor::getConnections()
Float_tgetCorrection(UInt_t flag = kQx)
Float_tgetCorrectionError(UInt_t flag = kQx)
Float_tgetCorrectionPhi(Float_t phi)
Float_tgetDirectivity()
virtual Option_t*TObject::GetDrawOption() const
static Long_tTObject::GetDtorOnly()
Float_tgetEt()
Float_tgetEventPlane(UInt_t EPcorr = kNoCorrection, UInt_t SubEvent = 0, UInt_t nHarmonic = 1)
TH1F*getEventPlaneCorrectionHist(UInt_t flag) const
Float_tgetEventPlaneParameter(UInt_t flag = kQx, Bool_t corr = kFALSE)
UInt_tgetEventPlaneParameterFromString(TString s)
Float_tgetEventPlaneWeight(UInt_t EPcorr = kNoCorrection, UInt_t SubEvent = 0, UInt_t nHarmonic = 1)
Float_tgetEventWeight()
Int_tgetFWCharge(HWallHitSim* wall)
vector<Int_t>getFWhits()
Int_tGetFWmoduleSize(HWallHitSim* wall)
Float_tgetFWSumChargeSpec()
Float_tgetFWSumZ(Int_t minZ = 1, Int_t maxZ = 99, UInt_t SubEvent = 0)
virtual const char*TObject::GetIconName() const
vector<TString>getLabelArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Int_tgetMdcWiresOuterMod()
virtual const char*TNamed::GetName() const
Int_tgetNbins(TString estimator)
Int_tgetNbins(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Int_tgetNumFWCells()
virtual char*TObject::GetObjectInfo(Int_t px, Int_t py) const
static Bool_tTObject::GetObjectStat()
virtual Option_t*TObject::GetOption() const
virtual HTask*HTask::getOwner()
Float_tgetRatioEtEz()
Float_tgetSelectedParticleCandCorrPerWire()
Float_tgetSelectedParticleCandSecNorm()
Float_tgetSmearValue(Int_t size = 1)
TStringgetStringCentralityClass(UInt_t e)
TStringgetStringCentralityEstimator(UInt_t e)
TStringgetStringEventPlaneCorrection(UInt_t e)
TStringgetStringEventPlaneParameter(UInt_t e)
TStringgetStringFWCutValues(UInt_t e)
TObjArrayHReconstructor::getTable()
virtual HTask*HReconstructor::getTask(const Char_t* name)
Float_tgetThetaWeight(HWallHitSim* wall, Float_t min = 3.5, Float_t max = 8.)
virtual const char*TNamed::GetTitle() const
virtual UInt_tTObject::GetUniqueID() const
Float_t*getUpEdgeArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Float_tgetVersion()
virtual Bool_tTObject::HandleTimer(TTimer* timer)
virtual ULong_tTNamed::Hash() const
virtual voidTObject::Info(const char* method, const char* msgfmt) const
virtual Bool_tTObject::InheritsFrom(const char* classname) const
virtual Bool_tTObject::InheritsFrom(const TClass* cl) const
virtual Bool_tinit()
virtual voidTObject::Inspect() constMENU
voidTObject::InvertBit(UInt_t f)
virtual TClass*IsA() const
virtual Bool_tTObject::IsEqual(const TObject* obj) const
virtual Bool_tHReconstructor::IsFolder() const
Bool_tisNewEvent()
Bool_tTObject::IsOnHeap() const
virtual Bool_tTNamed::IsSortable() const
virtual voidHReconstructor::isTimed(Bool_t flag = kTRUE)
Bool_tTObject::IsZombie() const
Bool_tloadCentralityEstimatorHist()
Int_tloadDayOfYear()
Bool_tloadEventPlaneCorrectionHist()
Bool_tloadFWCutValuesHist()
Bool_tloadParameterFile()
virtual voidTNamed::ls(Option_t* option = "") const
TH1F*makeClasses(TH1F* h, Float_t fractionXsection = 100., UInt_t centC = k10, Int_t direction = -1)
TH2D*makeEPresolution(TProfile2D* hist, Bool_t calcChi = kFALSE)
TH1*makeEPresolution(TH3* hist, Bool_t calcChi = kFALSE)
TH1F*makePercentiles(TH1F* hist, Float_t fractionXsection = 100., Int_t direction = -1)
voidTObject::MayNotUse(const char* method) const
Double_tModifiedBesselI(Int_t n, Double_t x) const
virtual HTask*HReconstructor::next(Int_t& errCode)
virtual Bool_tTObject::Notify()
voidTObject::Obsolete(const char* method, const char* asOfVers, const char* removedFromVers) const
static voidTObject::operator delete(void* ptr)
static voidTObject::operator delete(void* ptr, void* vp)
static voidTObject::operator delete[](void* ptr)
static voidTObject::operator delete[](void* ptr, void* vp)
void*TObject::operator new(size_t sz)
void*TObject::operator new(size_t sz, void* vp)
void*TObject::operator new[](size_t sz)
void*TObject::operator new[](size_t sz, void* vp)
HTask&HTask::operator=(const HTask&)
virtual voidTObject::Paint(Option_t* option = "")
Bool_tPassesCutsFW(HWallHitSim* wall)
Bool_tPassesCutsFW(HWallHitSim* wall, UInt_t eFWCut)
virtual voidTObject::Pop()
Int_tprint()
virtual voidTNamed::Print(Option_t* option = "") const
Bool_tprintCentralityClass(TString estimator)
Bool_tprintCentralityClass(TH1* htemp)
Bool_tprintCentralityClass(UInt_t centE = kTOFRPC, UInt_t centC = k10)
voidprintHitArray()
voidprintQVectors()
virtual voidHReconstructor::printTimer()
virtual Int_tTObject::Read(const char* name)
virtual voidTObject::RecursiveRemove(TObject* obj)
virtual Bool_treinit()
voidreset()
voidTObject::ResetBit(UInt_t f)
virtual voidHReconstructor::resetTimer()
virtual voidTObject::SaveAs(const char* filename = "", Option_t* option = "") constMENU
Bool_tsaveCentralityEstimatorHist()
Bool_tsaveEventPlaneCorrectionHist()
Bool_tsaveFWCutValuesHist()
Bool_tsaveParameterFile()
virtual voidTObject::SavePrimitive(ostream& out, Option_t* option = "")
voidHReconstructor::setActive(Bool_t state)MENU
voidTObject::SetBit(UInt_t f)
voidTObject::SetBit(UInt_t f, Bool_t set)
virtual Bool_tHTask::setConnections()
virtual voidTObject::SetDrawOption(Option_t* option = "")MENU
static voidTObject::SetDtorOnly(void* obj)
voidsetExcludeNoisyFWcells(Bool_t b = kTRUE)
voidHTask::setManual()
virtual voidTNamed::SetName(const char* name)MENU
virtual voidTNamed::SetNameTitle(const char* name, const char* title)
static voidTObject::SetObjectStat(Bool_t stat)
virtual voidHTask::setOwner(HTask* atask)
Bool_tsetParameterFile(TString ParameterFile)
voidsetPostFix(TString postFix)
voidsetReferenceMean(Float_t referenceMean)
virtual voidTNamed::SetTitle(const char* title = "")MENU
virtual voidTObject::SetUniqueID(UInt_t uid)
voidsetUseFWCut(UInt_t eFWCut, Bool_t b = kTRUE)
virtual voidShowMembers(TMemberInspector&)
virtual Int_tTNamed::Sizeof() const
virtual voidStreamer(TBuffer&)
voidStreamerNVirtual(TBuffer& ClassDef_StreamerNVirtual_b)
virtual voidTObject::SysError(const char* method, const char* msgfmt) const
Bool_tTObject::TestBit(UInt_t f) const
Int_tTObject::TestBits(UInt_t f) const
virtual voidTObject::UseCurrentStyle()
virtual voidTObject::Warning(const char* method, const char* msgfmt) const
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0)
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0) const
protected:
virtual voidTObject::DoError(int level, const char* location, const char* fmt, va_list va) const
voidTObject::MakeZombie()

Data Members

public:
enum eCentralityEstimator { kTOFRPC
kTOF
kRPC
kTOFRPCtot
kTOFtot
kRPCtot
kSelectedParticleCand
kSelectedParticleCandCorr
kSelectedParticleCandNorm
kSelectedParticleCandSecNorm
kSelectedParticleCandCorrPerWire
kPrimaryParticleCand
kMdcWires
kMdcWiresOuterMod
kFWSumChargeSpec
kFWSumChargeZ
kDirectivity
kRatioEtEz
kEt
kNumCentralityEstimator
};
enum eCentralityClass { k10
k5
k2
k1
k13
k15
k20
k1040
kFOPI
kNumCentralityClass
};
enum eEventPlaneCorrection { kNoCorrection
kSmear
kShiftFW
kScaleFW
kReCenter
kWidthEqualization
kFlatten
kWeightCharge
kWeightTheta
kScalingFW
kReCentering
kScaling
kFlattening
kRotation
kDefault
};
enum eEventPlaneParameter { kFWx
kFWy
kQx
kQy
kQxWCharge
kQyWCharge
kQx2
kQy2
kQx2WCharge
kQy2WCharge
kFourierC1
kFourierC2
kFourierC3
kFourierC4
kFourierC5
kFourierC6
kFourierC7
kFourierC8
kFourierS1
kFourierS2
kFourierS3
kFourierS4
kFourierS5
kFourierS6
kFourierS7
kFourierS8
kResolution
kResolutionWCharge
kChi
kNumEventPlaneParameter
};
enum eFWCut { kBetaCuts
kTimeCuts
kChargeCuts
kNumFWCutValues
};
enum TObject::EStatusBits { kCanDelete
kMustCleanup
kObjInCanvas
kIsReferenced
kHasUUID
kCannotPick
kNoContextMenu
kInvalidObject
};
enum TObject::[unnamed] { kIsOnHeap
kNotDeleted
kZombie
kBitMask
kSingleKey
kOverwrite
kWriteDelete
};
protected:
Bool_tHReconstructor::fActiveActive flag
TList*HReconstructor::fHistogramsList of histograms generated by this reconstructor.
Bool_tHTask::fIsTimedtimer flag
TStringTNamed::fNameobject identifier
TObjArrayHReconstructor::fOutputs
TStopwatchHReconstructor::fTimerTask timer
TStringTNamed::fTitleobject title
Bool_tHTask::isConnected
Bool_tHTask::isInitialised
Bool_tHTask::manual
HTask*HTask::owner
private:
vector<SimpleQVector*>arrayOfHits
UInt_tcurrentEventSeqNumberrunning seq. num. of event
Bool_texcludeNoisyFWcells
HCategory*fCatWallHitForwardWall category
vector<vector<TH1*> >fCentralityHist
vector<TH1*>fCentralityPercentileHist
UInt_tfDayOfYearday of the year extracted fron filename (sim=0)
vector<TH1*>fEstimatorHist
UInt_tfEventPlaneCorrectionFlag
vector<TH2*>fEventPlaneCorrectionHist
vector<Float_t>fFWChargeCuts
vector<vector<TH1*> >fFWCutValuesHist
vector<Float_t>fFWmaxBeta
vector<Float_t>fFWminBeta
vector<Float_t>fFWminCharge
TFile*fFile
TStringfParameterFilelocation and name of ParameterFile
HCategory*fParticleCandCatParticleCand category
HCategory*fParticleEvtInfoCatParticleEvtInfo category
TStringfPostFix
Bool_tfQVectorCalcDone
TRandom*fRandom
Float_tfReferenceMeanSelTracktrack running mean reference point
vector<Int_t>iFWHitvectorlist of selected and shuffled FW hit-num.
Bool_tisSimulationfor simulation if catGeantKine is available
vector<Bool_t>useFWCut
vector<Float_t>vQPhi

Class Charts

Inheritance Inherited Members Includes Libraries
Class Charts

Function documentation

HParticleEvtChara(const Text_t* name = "HParticleEvtChara", const Text_t* title = "HParticleEvtChara")
~HParticleEvtChara()
Bool_t init()
Bool_t reinit()
Int_t loadDayOfYear()
Bool_t setParameterFile(TString ParameterFile)
Bool_t loadParameterFile()
 read parameter file
Bool_t saveParameterFile()
 read parameter file
Bool_t loadCentralityEstimatorHist()
Bool_t saveCentralityEstimatorHist()
Bool_t loadEventPlaneCorrectionHist()
Bool_t addEventPlaneCorrectionHist(TProfile2D* hist, UInt_t epParam = kQx)
Bool_t saveEventPlaneCorrectionHist()
Bool_t loadFWCutValuesHist()
Bool_t addFWCutValuesHist(TH1* hist, Int_t cell, UInt_t eFWCut)
Bool_t saveFWCutValuesHist()
Int_t execute()
Int_t print()
void reset()
 Reset.
Bool_t isNewEvent()
Float_t getEventWeight()
Int_t GetFWmoduleSize(HWallHitSim* wall)
Bool_t PassesCutsFW(HWallHitSim* wall)
Bool_t PassesCutsFW(HWallHitSim* wall, UInt_t eFWCut)
Int_t getFWCharge(HWallHitSim* wall)
TH1F* getEventPlaneCorrectionHist(UInt_t flag) const
Float_t getCorrection(UInt_t flag = kQx)
Float_t getCorrectionError(UInt_t flag = kQx)
Float_t getSmearValue(Int_t size = 1)
Float_t getThetaWeight(HWallHitSim* wall, Float_t min = 3.5, Float_t max = 8.)
void fillHitArray()
fill Hits into HitArray for QVector
void printHitArray()
vector<Int_t> getFWhits()
Bool_t fillQVectors(UInt_t EPcorr = kNoCorrection, UInt_t nHarmonic = 1)
HWallHit  into QVector
void printQVectors()
Float_t getCorrectionPhi(Float_t phi)
used Formular in Footnoot1 of Phys.Rev. C58 (1998) 1671-1678
also appendix of Phys.Rev. C56 (1997) 3254-3264
Float_t getEventPlane(UInt_t EPcorr = kNoCorrection, UInt_t SubEvent = 0, UInt_t nHarmonic = 1)
Float_t getEventPlaneWeight(UInt_t EPcorr = kNoCorrection, UInt_t SubEvent = 0, UInt_t nHarmonic = 1)
Float_t getEventPlaneParameter(UInt_t flag = kQx, Bool_t corr = kFALSE)
Int_t getCentralityClass(TString estimator)
 legacy code
 Return centrality class, default in 5% of total cross section with estimator
 or with preset classes like FOPI {6.3%, 21.0%, 30.9%}
Float_t getCentralityEstimator(UInt_t centE = kTOFRPC)
TH1F* getCentralityClassHist(UInt_t centE = kTOFRPC, UInt_t centC = k10) const
TH1F* getCentralityPercentileHist(UInt_t centE = kTOFRPC) const
Int_t getCentralityClass(UInt_t centE = kTOFRPC, UInt_t centC = k10)
 Return centrality class, default in 5% of total cross section with estimator
 or with preset classes like FOPI {6.3%, 21.0%, 30.9%}
Bool_t printCentralityClass(TString estimator)
 legacy code
 print all CentralityClasses in the Estimator
Bool_t printCentralityClass(UInt_t centE = kTOFRPC, UInt_t centC = k10)
 print all CentralityClasses in the Estimator
Bool_t printCentralityClass(TH1* htemp)
 print all CentralityClasses in the Estimator
TH2D* makeEPresolution(TProfile2D* hist, Bool_t calcChi = kFALSE)
TH1* makeEPresolution(TH3* hist, Bool_t calcChi = kFALSE)
Double_t ModifiedBesselI(Int_t n, Double_t x) const
 compute half-integer modified Bessel functions
 order: n>0, for n>5/2, interpolation is used (improve this by using recurrence!!!)
Double_t ComputeResolution(Double_t x, Int_t n = 1) const
 Computes resolution for Event Plane method
Double_t FindXi(Double_t res, Double_t prec = 1e-6, Int_t n = 1) const
 Computes x(res) for Event Plane method
Bool_t addEstimatorHist(TH1F* hist, Float_t fractionXsection = 100., UInt_t centE = kTOFRPC, Int_t direction = -1)
TH1F* makePercentiles(TH1F* hist, Float_t fractionXsection = 100., Int_t direction = -1)
TH1F* makeClasses(TH1F* h, Float_t fractionXsection = 100., UInt_t centC = k10, Int_t direction = -1)
Int_t getNbins(TString estimator)
 legacy code
 Number of Bins (CentralityClasses and Over- and Underflow) in the Estimator
Int_t getNbins(UInt_t centE = kTOFRPC, UInt_t centC = k10)
 Number of Bins (CentralityClasses and Over- and Underflow) in the Estimator
Int_t getCentralityClassNbins(UInt_t centC = k10)
Float_t* getCentralityClassArray(UInt_t centC = k10)
Float_t* getUpEdgeArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Float_t* getBinCenterArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
vector<TString> getLabelArray(UInt_t centE = kTOFRPC, UInt_t centC = k10)
Int_t getCentralityClass5(TString estimator)
 Return centrality class 5% of total cross section
Int_t getCentralityClass10(TString estimator)
 Return centrality class 5% of total cross section
Float_t getCentralityPercentile(TString estimator)
 legacy code
Float_t getCentralityPercentile(UInt_t centE = kTOFRPC)
Int_t getMdcWiresOuterMod()
Float_t getSelectedParticleCandSecNorm()
Float_t getSelectedParticleCandCorrPerWire()
Float_t getFWSumChargeSpec()
Float_t getFWSumZ(Int_t minZ = 1, Int_t maxZ = 99, UInt_t SubEvent = 0)
Float_t getEt()
Float_t getRatioEtEz()
Float_t getDirectivity()
Int_t getCentralityEstimatorBinSize(UInt_t e)
Int_t getCentralityEstimatorBinMax(UInt_t e)
TString getStringCentralityEstimator(UInt_t e)
UInt_t getCentralityEstimatorFromString(TString s)
TString getStringCentralityClass(UInt_t e)
Float_t getCentralityClassBinSize(UInt_t e)
TString getStringEventPlaneCorrection(UInt_t e)
TString getStringEventPlaneParameter(UInt_t e)
UInt_t getEventPlaneParameterFromString(TString s)
TString getStringFWCutValues(UInt_t e)
Float_t getVersion()
{return 0.7;}
Int_t getNumFWCells()
{return MAXFWCELLS;}
Bool_t finalize()
{ return kTRUE;}
void setPostFix(TString postFix)
setter:
{fPostFix = postFix;}
void setUseFWCut(UInt_t eFWCut, Bool_t b = kTRUE)
{useFWCut[eFWCut] = b;}
void setExcludeNoisyFWcells(Bool_t b = kTRUE)
void setReferenceMean(Float_t referenceMean)
{fReferenceMeanSelTrack = referenceMean;}
Bool_t isFlagSet(UInt_t flag, UInt_t status)
{ return (flag==(status&flag));}
Float_t getPhi(Float_t x, Float_t y)
{ return TVector2::Phi_0_2pi( TMath::Pi()+TMath::ATan2(-y,-x) );}