Commit 351e3d68 authored by Semen Lebedev's avatar Semen Lebedev
Browse files

Major update in LMVM

parent 78ee7a1e
......@@ -47,24 +47,21 @@ ${Boost_LIBRARY_DIRS}
link_directories( ${LINK_DIRECTORIES})
set(SRCS
CbmAnaDielectronTask.cxx
CbmAnaDielectronTaskDraw.cxx
CbmAnaDielectronTaskDrawAll.cxx
CbmAnaLmvmDrawStudy.cxx
CbmLmvmHist.cxx
CbmAnaDielectronStudyReportAll.cxx
CbmAnaDielectronReports.cxx
CbmHaddBase.cxx
LmvmTask.cxx
LmvmDraw.cxx
LmvmDrawAll.cxx
LmvmHist.cxx
LmvmUtils.cxx
)
IF (SSE_FOUND)
Message(STATUS "Analysis will be compiled with SSE support")
Message(STATUS "LMVM analysis will be compiled with SSE support")
ADD_DEFINITIONS(-DHAVE_SSE)
SET_SOURCE_FILES_PROPERTIES(${SRCS} PROPERTIES COMPILE_FLAGS
"-msse -O3")
ELSE (SSE_FOUND)
MESSAGE(STATUS "Analysis will be compiled without SSE support")
MESSAGE(STATUS "LMVM analysis will be compiled without SSE support")
SET_SOURCE_FILES_PROPERTIES(${SRCS} PROPERTIES COMPILE_FLAGS
"-O3")
ENDIF (SSE_FOUND)
......@@ -72,10 +69,10 @@ ENDIF (SSE_FOUND)
ADD_DEFINITIONS(-DDO_TPCCATRACKER_EFF_PERFORMANCE -DNonhomogeneousField -DCBM -DUSE_TIMERS)
Set(LINKDEF CbmDiElectronAnalysisLinkDef.h)
Set(LIBRARY_NAME CbmDiElectronAnalysis)
Set(LINKDEF LmvmLinkDef.h)
Set(LIBRARY_NAME CbmLmvm)
Set(DEPENDENCIES
Littrack KF L1 CbmRichBase CbmRichReco CbmRecoBase CbmBase CbmData)
Littrack LittrackQA KF L1 CbmRichBase CbmRichReco CbmRecoBase CbmBase CbmData)
Set(DEFINITIONS -DDO_TPCCATRACKER_EFF_PERFORMANCE -DNonhomogeneousField -DCBM -DUSE_TIMERS)
GENERATE_LIBRARY()
......
/* Copyright (C) 2015 Justus-Liebig-Universitaet Giessen, Giessen
/* Copyright (C) 2015-2021 Justus-Liebig-Universitaet Giessen, Giessen
SPDX-License-Identifier: GPL-3.0-only
Authors: Elena Lebedeva [committer] */
Authors: Elena Lebedeva [committer], Semen Lebedev */
/**
* @author Elena Lebedeva <e.lebedeva@gsi.de>
* @since 2015
* @version 1.0
**/
#ifndef CBM_LMVM_KINE_PARAMS_H
#define CBM_LMVM_KINE_PARAMS_H
#ifndef CBM_LMVM_KINEMATIC_PARAMS_H
#define CBM_LMVM_KINEMATIC_PARAMS_H
#include "CbmLmvmCandidate.h"
#include "CbmMCTrack.h"
#include "TLorentzVector.h"
#include "TMath.h"
#include "LmvmCand.h"
#define M2E 2.6112004954086e-7
class CbmLmvmKinematicParams {
public:
Double_t fMomentumMag; // Absolute value of momentum
Double_t fPt; // Transverse momentum
Double_t fRapidity; // Rapidity
Double_t fMinv; // Invariant mass
Double_t fAngle; // Opening angle
Double_t fMomentumMag = 0.; // Absolute value of momentum
Double_t fPt = 0.; // Transverse momentum
Double_t fRapidity = 0.; // Rapidity
Double_t fMinv = 0.; // Invariant mass
Double_t fAngle = 0.; // Opening angle
/*
* Calculate kinematic parameters for MC tracks.
*/
static CbmLmvmKinematicParams KinematicParamsWithMcTracks(const CbmMCTrack* mctrackP, const CbmMCTrack* mctrackM)
static CbmLmvmKinematicParams Create(const CbmMCTrack* mctrackP, const CbmMCTrack* mctrackM)
{
CbmLmvmKinematicParams params;
if (mctrackP == nullptr || mctrackM == nullptr) return params;
TVector3 momP; //momentum e+
mctrackP->GetMomentum(momP);
......@@ -64,10 +60,10 @@ public:
/*
* Calculate kinematic parameters for LMVM candidates.
*/
static CbmLmvmKinematicParams KinematicParamsWithCandidates(const CbmLmvmCandidate* candP,
const CbmLmvmCandidate* candM)
static CbmLmvmKinematicParams Create(const LmvmCand* candP, const LmvmCand* candM)
{
CbmLmvmKinematicParams params;
if (candP == nullptr || candM == nullptr) return params;
TLorentzVector lorVecP(candP->fMomentum, candP->fEnergy);
TLorentzVector lorVecM(candM->fMomentum, candM->fEnergy);
......
/* Copyright (C) 2015 Justus-Liebig-Universitaet Giessen, Giessen
/* Copyright (C) 2015-2021 Justus-Liebig-Universitaet Giessen, Giessen
SPDX-License-Identifier: GPL-3.0-only
Authors: Elena Lebedeva [committer] */
Authors: Elena Lebedeva [committer], Semen Lebedev */
/**
* @author Elena Lebedeva <e.lebedeva@gsi.de>
* @since 2015
* @version .0
**/
#ifndef CBM_LMVM_CANDIDATE_H
#define CBM_LMVM_CANDIDATE_H
#ifndef LMVM_CAND_H
#define LMVM_CAND_H
#include "TVector3.h"
class CbmLmvmCandidate {
#include "LmvmDef.h"
class LmvmCand {
public:
CbmLmvmCandidate()
: fPosition()
, fMomentum()
, fMass(0.)
, fEnergy(0.)
, fRapidity(0.)
, fCharge(0)
, fChi2Prim(0.)
, fChi2sts(0.)
, fMcMotherId(-1)
, fEventNumber(0.)
, fStsMcTrackId(-1)
, fRichMcTrackId(-1)
, fTrdMcTrackId(-1)
, fTofMcTrackId(-1)
, fStsInd(-1)
, fRichInd(-1)
, fTrdInd(-1)
, fTofInd(-1)
, fIsElectron(kFALSE)
, fIsMcSignalElectron(kFALSE)
, fIsMcPi0Electron(kFALSE)
, fIsMcGammaElectron(kFALSE)
, fIsMcEtaElectron(kFALSE)
, fMcPdg(-1)
, fIsGamma(kFALSE)
, fDSts(0.)
, fIsTtCutElectron(kFALSE)
, fIsStCutElectron(kFALSE)
, fIsRtCutElectron(kFALSE)
, fIsMvd1CutElectron(kFALSE)
, fIsMvd2CutElectron(kFALSE)
, fRichAnn(0.)
, fTrdAnn(0.)
, fMass2(0.)
LmvmCand() {}
void ResetMcParams()
{
fMcSrc = ELmvmSrc::Undefined;
fMcMotherId = -1;
fStsMcTrackId = -1;
fRichMcTrackId = -1;
fTrdMcTrackId = -1;
fTofMcTrackId = -1;
}
void ResetMcParams()
void SetIsTopologyCutElectron(ELmvmTopologyCut cut, bool value)
{
fIsMcSignalElectron = false;
fIsMcPi0Electron = false;
fIsMcGammaElectron = false;
fIsMcEtaElectron = false;
fMcMotherId = -1;
fStsMcTrackId = -1;
fRichMcTrackId = -1;
fTrdMcTrackId = -1;
fTofMcTrackId = -1;
if (cut == ELmvmTopologyCut::TT) { fIsTtCut = value; }
if (cut == ELmvmTopologyCut::ST) { fIsStCut = value; }
if (cut == ELmvmTopologyCut::RT) { fIsRtCut = value; }
}
bool IsCutTill(ELmvmAnaStep step) const
{
if (step == ELmvmAnaStep::Mc || step == ELmvmAnaStep::Acc || step == ELmvmAnaStep::Reco) return true;
if (step == ELmvmAnaStep::Chi2Prim && fIsChi2Prim) return true;
if (step == ELmvmAnaStep::ElId && IsCutTill(ELmvmAnaStep::Chi2Prim) && fIsElectron) return true;
if (step == ELmvmAnaStep::GammaCut && IsCutTill(ELmvmAnaStep::ElId) && fIsGammaCut) return true;
if (step == ELmvmAnaStep::Mvd1Cut && IsCutTill(ELmvmAnaStep::GammaCut) && fIsMvd1Cut) return true;
if (step == ELmvmAnaStep::Mvd2Cut && IsCutTill(ELmvmAnaStep::Mvd1Cut) && fIsMvd2Cut) return true;
if (step == ELmvmAnaStep::StCut && IsCutTill(ELmvmAnaStep::Mvd2Cut) && fIsStCut) return true;
if (step == ELmvmAnaStep::RtCut && IsCutTill(ELmvmAnaStep::StCut) && fIsRtCut) return true;
if (step == ELmvmAnaStep::TtCut && IsCutTill(ELmvmAnaStep::RtCut) && fIsTtCut) return true;
if (step == ELmvmAnaStep::PtCut && IsCutTill(ELmvmAnaStep::TtCut) && fIsPtCut) return true;
return false;
}
// track parameters
TVector3 fPosition;
TVector3 fMomentum;
Double_t fMass;
Double_t fEnergy;
Double_t fRapidity;
Int_t fCharge;
Double_t fChi2Prim;
Double_t fChi2sts;
double fMass = 0.;
double fEnergy = 0.;
double fRapidity = 0.;
int fCharge = 0;
double fChi2Prim = 0.;
double fChi2sts = 0.;
int fMcMotherId = -1;
int fEventNumber = 0;
int fStsMcTrackId = -1;
int fRichMcTrackId = -1;
int fTrdMcTrackId = -1;
int fTofMcTrackId = -1;
int fStsInd = -1;
int fRichInd = -1;
int fTrdInd = -1;
int fTofInd = -1;
int fMcPdg = -1;
double fRichAnn = 0.;
double fTrdAnn = 0.;
double fMass2 = 0.;
Int_t fMcMotherId;
Int_t fEventNumber;
Int_t fStsMcTrackId;
Int_t fRichMcTrackId;
Int_t fTrdMcTrackId;
Int_t fTofMcTrackId;
Int_t fStsInd;
Int_t fRichInd;
Int_t fTrdInd;
Int_t fTofInd;
Bool_t fIsElectron;
Bool_t fIsMcSignalElectron;
Bool_t fIsMcPi0Electron;
Bool_t fIsMcGammaElectron;
Bool_t fIsMcEtaElectron;
// Cuts. If true then cut is passed
bool fIsChi2Prim = false;
bool fIsElectron = false;
bool fIsGammaCut = false;
bool fIsMvd1Cut = false;
bool fIsMvd2Cut = false;
bool fIsTtCut = false;
bool fIsStCut = false;
bool fIsRtCut = false;
bool fIsPtCut = false;
Int_t fMcPdg;
Bool_t fIsGamma;
Double_t fDSts;
Bool_t fIsTtCutElectron;
Bool_t fIsStCutElectron;
Bool_t fIsRtCutElectron;
Bool_t fIsMvd1CutElectron;
Bool_t fIsMvd2CutElectron;
Double_t fRichAnn;
Double_t fTrdAnn;
Double_t fMass2;
// MC
ELmvmSrc fMcSrc = ELmvmSrc::Undefined;
bool IsMcSignal() const { return fMcSrc == ELmvmSrc::Signal; }
bool IsMcPi0() const { return fMcSrc == ELmvmSrc::Pi0; }
bool IsMcGamma() const { return fMcSrc == ELmvmSrc::Gamma; }
bool IsMcEta() const { return fMcSrc == ELmvmSrc::Eta; }
// for candidates BG is all candidates which are not signal
bool fIsMcBg() const { return fMcSrc != ELmvmSrc::Signal; }
};
#endif
/* Copyright (C) 2015-2018 Justus-Liebig-Universitaet Giessen, Giessen
/* Copyright (C) 2015-2021 Justus-Liebig-Universitaet Giessen, Giessen
SPDX-License-Identifier: GPL-3.0-only
Authors: Elena Lebedeva [committer], Gregor Pitsch */
Authors: Elena Lebedeva [committer], Gregor Pitsch, Semen Lebedev */
/**
* @author Elena Lebedeva <e.lebedeva@gsi.de>
* @since 2015
* @version 1.0
**/
#ifndef LMVM_CUTS_H
#define LMVM_CUTS_H
#include <Logger.h>
#ifndef CBM_LMVM_CUTS_H
#define CBM_LMVM_CUTS_H
#include <math.h>
#include "TObject.h"
#include "LmvmDef.h"
#include <iostream>
class CbmLmvmCuts {
class LmvmCuts {
public:
CbmLmvmCuts()
: fMomentumCut(0.)
, fChiPrimCut(0.)
, fPtCut(0.)
, fAngleCut(0.)
, fGammaCut(0.)
, fStCutAngle(0.)
, fStCutPP(0.)
, fTtCutAngle(0.)
, fTtCutPP(0.)
, fRtCutAngle(0.)
, fRtCutPP(0.)
, fMvd1CutP(0.)
, fMvd1CutD(0.)
, fMvd2CutP(0.)
, fMvd2CutD(0.)
LmvmCuts() {}
bool IsTopologyCutOk(ELmvmTopologyCut cut, double mom1, double mom2, double minAngle)
{
SetDefaultCuts();
double angleCut = 0., ppCut = 0.;
if (cut == ELmvmTopologyCut::ST) {
angleCut = fStCutAngle;
ppCut = fStCutPP;
}
else if (cut == ELmvmTopologyCut::RT) {
angleCut = fRtCutAngle;
ppCut = fRtCutPP;
}
else if (cut == ELmvmTopologyCut::TT) {
angleCut = fTtCutAngle;
ppCut = fTtCutPP;
}
else {
LOG(error) << "LmvmCuts::IsTopologyCut cut is not defined.";
}
double sqrt_mom = std::sqrt(mom1 * mom2);
double val = -1. * (angleCut / ppCut) * sqrt_mom + angleCut;
if (!(sqrt_mom < ppCut && val > minAngle)) return true;
return false;
}
/*
* Set default electron ID and analysis cuts.
*/
void SetDefaultCuts()
{
//electron ID cuts, we use CbmLitGlobalElectronId for identification
fMomentumCut = -1.; // if cut < 0 it is not used
bool IsChi2PrimaryOk(double chi2Prim) { return (chi2Prim < fChi2PrimCut); }
// analysis cuts auau
fPtCut = 0.2;
fAngleCut = 1.;
fChiPrimCut = 3.;
fGammaCut = 0.025;
//fStCutAngle = 1.5;
//fStCutPP = 1.5;
fStCutAngle = 2.4;
fStCutPP = 1.;
//fTtCutAngle = 0.75;
//fTtCutPP = 4.0;
fTtCutAngle = 2.0;
fTtCutPP = 3.0;
//fRtCutAngle = 1.0;
//fRtCutPP = 2.5;
fRtCutAngle = 1.2;
fRtCutPP = 1.6;
fMvd1CutP = 1.2;
fMvd1CutD = 0.4;
fMvd2CutP = 1.5;
fMvd2CutD = 0.5;
bool IsGammaCutOk(double minv) { return (minv >= fGammaCut); }
// analysis cuts agag
/* fPtCut = 0.2;
fAngleCut = 1.;
fChiPrimCut = 3.;
fGammaCut = 0.025;
fStCutAngle = 2.4;
fStCutPP = 1.;
fTtCutAngle = 1.5;
fTtCutPP = 1.7;
fRtCutAngle = 1.2;
fRtCutPP = 1.6;
fMvd1CutP = 1.2;
fMvd1CutD = 0.4;
fMvd2CutP = 1.5;
fMvd2CutD = 0.5;
*/
bool IsPtCutOk(double pt) { return (pt > fPtCut); }
bool IsMvdCutOk(int stationNum, double dmvd, double mom)
{
if (stationNum <= 0 || stationNum > 2) {
LOG(error) << "LmvmCuts::IsMvdCut stationNum is not in valid. stationNum = " << stationNum;
return false;
}
// it is assumed that stationNum can be 1 or 2
double cutD = (stationNum == 1) ? fMvd1CutD : fMvd2CutD;
double cutP = (stationNum == 1) ? fMvd1CutP : fMvd2CutP;
double val = -1. * (cutP / cutD) * dmvd + cutP;
if (!(dmvd < cutD && val > mom)) return true;
return false;
}
/*
* Print out cuts.
*/
void Print()
std::string ToString()
{
std::cout << "Used cuts:" << std::endl
<< "fChiPrimCut = " << fChiPrimCut << std::endl
<< "fPtCut = " << fPtCut << std::endl
<< "fAngleCut = " << fAngleCut << std::endl
<< "fGammaCut = " << fGammaCut << std::endl
<< "fStCut (ang,pp) = (" << fStCutAngle << "," << fStCutPP << ")" << std::endl
<< "fRtCut (ang,pp) = (" << fRtCutAngle << "," << fRtCutPP << ")" << std::endl
<< "fTtCut (ang,pp) = (" << fTtCutAngle << "," << fTtCutPP << ")" << std::endl
<< "fMvd1Cut (p,d) = (" << fMvd1CutP << "," << fMvd1CutD << ")" << std::endl
<< "fMvd2Cut (p,d) = (" << fMvd2CutP << "," << fMvd2CutD << ")" << std::endl
<< "fMomentumCut = " << fMomentumCut << std::endl;
std::stringstream ss;
ss << "LMVM cuts:" << std::endl
<< "fChiPrimCut = " << fChi2PrimCut << std::endl
<< "fPtCut = " << fPtCut << std::endl
<< "fAngleCut = " << fAngleCut << std::endl
<< "fGammaCut = " << fGammaCut << std::endl
<< "fStCut (ang,pp) = (" << fStCutAngle << "," << fStCutPP << ")" << std::endl
<< "fRtCut (ang,pp) = (" << fRtCutAngle << "," << fRtCutPP << ")" << std::endl
<< "fTtCut (ang,pp) = (" << fTtCutAngle << "," << fTtCutPP << ")" << std::endl
<< "fMvd1Cut (p,d) = (" << fMvd1CutP << "," << fMvd1CutD << ")" << std::endl
<< "fMvd2Cut (p,d) = (" << fMvd2CutP << "," << fMvd2CutD << ")" << std::endl
<< "fMomentumCut = " << fMomentumCut << std::endl;
return ss.str();
}
public:
// ID cuts, we use CbmLitGlobalElectronId for identification
Double_t fMomentumCut; // if cut < 0 then it will not be used
double fMomentumCut = -1.; // if fMomentumCut < 0 it is not used
// Analysis cuts
Double_t fChiPrimCut;
Double_t fPtCut;
Double_t fAngleCut;
Double_t fGammaCut;
Double_t fStCutAngle;
Double_t fStCutPP;
Double_t fTtCutAngle;
Double_t fTtCutPP;
Double_t fRtCutAngle;
Double_t fRtCutPP;
Double_t fMvd1CutP;
Double_t fMvd1CutD;
Double_t fMvd2CutP;
Double_t fMvd2CutD;
double fPtCut = 0.2;
double fAngleCut = 1.;
double fChi2PrimCut = 3.;
double fGammaCut = 0.025;
double fStCutAngle = 2.4;
double fStCutPP = 1.;
double fTtCutAngle = 2.0;
double fTtCutPP = 3.0;
double fRtCutAngle = 1.2;
double fRtCutPP = 1.6;
double fMvd1CutP = 1.2;
double fMvd1CutD = 0.4;
double fMvd2CutP = 1.5;
double fMvd2CutD = 0.5;
};
#endif
/* Copyright (C) 2021 Justus-Liebig-Universitaet Giessen, Giessen
SPDX-License-Identifier: GPL-3.0-only
Authors: Semen Lebedev [committer] */
#ifndef LMVM_DEF_H
#define LMVM_DEF_H
#include "Rtypes.h"
#include <string>
class TH1D;
enum class ELmvmTopologyCut : int
{
ST,
RT,
TT
};
enum class ELmvmSrc : int
{
Signal = 0,
Bg = 1,
Pi0 = 2,
Gamma = 3,
Eta = 4,
Undefined = 5
};
enum class ELmvmAnaStep : int
{
Mc = 0,
Acc = 1,
Reco = 2,
Chi2Prim = 3,
ElId = 4,
GammaCut = 5,
Mvd1Cut = 6,
Mvd2Cut = 7,
StCut = 8,
RtCut = 9,
TtCut = 10,
PtCut = 11,
Undefined = 12
};
enum class ELmvmBgPairSrc : int
{
GG = 0, // gamma-gamma
PP = 1, // pi0-pi0
OO = 2, // other-other
GP = 3, // gamma-pi0
GO = 4, // gamma-other
PO = 5, // pi0-other
Undefined = 6
};
enum class ELmvmSignal : int
{
Inmed = 0,
Qgp = 1,
Omega = 2,
Phi = 3,
OmegaD = 4
};
class LmvmDataXYInd {
public:
LmvmDataXYInd() {}
LmvmDataXYInd(double x, double y, int ind) : fX(x), fY(y), fInd(ind) {}
double fX = 0.;
double fY = 0.;
int fInd = 0;
};
class LmvmDataAngMomInd {
public:
LmvmDataAngMomInd() {}
LmvmDataAngMomInd(double angle, double mom, int ind) : fAngle(angle), fMom(mom), fInd(ind) {}
double fAngle = 0.;
double fMom = 0.;
int fInd = 0;
};
class LmvmDrawMinvData {
public:
LmvmDrawMinvData() {}
LmvmDrawMinvData(TH1D* h, Color_t fillColor, Color_t lineColor, int lineWidth, Style_t fillStyle,
const std::string& legend)
: fH(h)
, fFillColor(fillColor)
, fLineColor(lineColor)
, fLineWidth(lineWidth)
, fFillStyle(fillStyle)