Skip to content

Commit 6288a38

Browse files
sofiatomassinisofiatomassini
committed
CF task
1 parent ff43312 commit 6288a38

File tree

5 files changed

+649
-332
lines changed

5 files changed

+649
-332
lines changed

‎PWGCF/DataModel/femtotest.h‎

Lines changed: 187 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,187 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
//
12+
/// \brief create a table for single track selection.
13+
/// \author Sofia Tomassini
14+
/// \since 30 May 2023
15+
16+
#ifndef PWGCF_DATAMODEL_FEMTOTEST_H_
17+
#define PWGCF_DATAMODEL_FEMTOTEST_H_
18+
19+
//#include "Framework/ASoA.h"
20+
//#include "Framework/DataTypes.h"
21+
//#include "Framework/AnalysisDataModel.h"
22+
//#include "Common/DataModel/PIDResponse.h"
23+
//#include "Framework/Logger.h"
24+
//#include "Common/DataModel/Multiplicity.h"
25+
26+
#include <vector>
27+
#include "TLorentzVector.h"
28+
#include "TVector3.h"
29+
30+
void comb(std::vector<std::vector<int>> &indxs, int N, int K)
31+
{
32+
std::string bitmask(K, 1); // K leading 1's
33+
bitmask.resize(N, 0); // N-K trailing 0's
34+
35+
// print integers and permute bitmask
36+
do {
37+
std::vector<int> temp;
38+
for (int i = 0; i < N; ++i) // [0..N-1] integers
39+
{
40+
if (bitmask[i]){
41+
temp.push_back(i);
42+
}
43+
}
44+
indxs.push_back(temp);
45+
} while (std::prev_permutation(bitmask.begin(), bitmask.end()));
46+
}
47+
48+
//====================================================================================
49+
50+
class FemtoParticle {
51+
public:
52+
FemtoParticle(){}
53+
FemtoParticle(TLorentzVector* fourmomentum, const float& eta, const float& phi);
54+
//FemtoParticle(TLorentzVector* fourmomentum);
55+
FemtoParticle(const float& E, const float& px, const float& py, const float& pz, const float& eta, const float& phi);
56+
//FemtoParticle(const float& E, const float& px, const float& py, const float& pz);
57+
FemtoParticle(const FemtoParticle& obj);
58+
FemtoParticle(const FemtoParticle* obj);
59+
~FemtoParticle();
60+
FemtoParticle& operator=(const FemtoParticle &obj);
61+
62+
void Set4momentum(TLorentzVector* fourmomentum){_fourmomentum = fourmomentum;}
63+
void SetEta(const float& eta){_eta = eta;}
64+
void SetPhi(const float& phi){_phi = phi;}
65+
void SetSign(const float& sign){_sign = sign;}
66+
void SetMagField(const float& magField){_magField = magField;}
67+
68+
TLorentzVector* Get4momentum() const {return _fourmomentum;}
69+
float GetEta() const {return _eta;}
70+
float GetPhi() const {return _phi;}
71+
float GetPhiStar(const float& radius = 1.2){return _phi + asin( -0.3*_magField*_sign*radius/(2.0*_fourmomentum->Pt()) );}
72+
float GetMagField() const {return _magField;}
73+
74+
private:
75+
float _eta, _phi, _sign, _magField;
76+
TLorentzVector* _fourmomentum;
77+
};
78+
79+
FemtoParticle::FemtoParticle(TLorentzVector* fourmomentum, const float& eta, const float& phi){
80+
_fourmomentum = fourmomentum;
81+
_eta = eta;
82+
_phi = phi;
83+
}
84+
85+
//FemtoParticle::FemtoParticle(TLorentzVector* fourmomentum){
86+
// _fourmomentum = fourmomentum;
87+
//}
88+
89+
FemtoParticle::FemtoParticle(const float& E, const float& px, const float& py, const float& pz, const float& eta, const float& phi){
90+
_fourmomentum = new TLorentzVector(px, py, pz, E);
91+
_eta = eta;
92+
_phi = phi;
93+
}
94+
95+
//FemtoParticle::FemtoParticle(const float& E, const float& px, const float& py, const float& pz){
96+
// _fourmomentum = new TLorentzVector(*px, *py, *pz, *E);
97+
//}
98+
99+
FemtoParticle::FemtoParticle(const FemtoParticle& obj){
100+
Set4momentum(obj.Get4momentum());
101+
SetEta(obj.GetEta());
102+
SetPhi(obj.GetPhi());
103+
}
104+
105+
FemtoParticle::FemtoParticle(const FemtoParticle* obj){
106+
Set4momentum(obj->Get4momentum());
107+
SetEta(obj->GetEta());
108+
SetPhi(obj->GetPhi());
109+
}
110+
111+
FemtoParticle::~FemtoParticle()
112+
{
113+
}
114+
115+
FemtoParticle& FemtoParticle::operator=(const FemtoParticle &obj)
116+
{
117+
if (this != &obj) {
118+
Set4momentum(obj.Get4momentum());
119+
SetEta(obj.GetEta());
120+
SetPhi(obj.GetPhi());
121+
}
122+
123+
return *this;
124+
}
125+
126+
//====================================================================================
127+
128+
129+
class FemtoPair{
130+
public:
131+
FemtoPair(){};
132+
FemtoPair(FemtoParticle* first, FemtoParticle* second){_first = first; _second = second;}
133+
FemtoPair(FemtoParticle* first, FemtoParticle* second, const bool& isidentical){_first = first; _second = second; _isidentical = isidentical;}
134+
135+
FemtoPair(const FemtoPair& obj){ SetFirstParticle(obj.GetFirstParticle()); SetSecondParticle(obj.GetSecondParticle()); }
136+
FemtoPair(const FemtoPair* obj){ SetFirstParticle(obj->GetFirstParticle()); SetSecondParticle(obj->GetSecondParticle()); }
137+
~FemtoPair(){}
138+
FemtoPair& operator=(const FemtoPair &obj){ if (this != &obj) {SetFirstParticle(obj.GetFirstParticle()); SetSecondParticle(obj.GetSecondParticle());}
139+
return *this;
140+
}
141+
142+
void SetFirstParticle(FemtoParticle* first){_first = first;}
143+
void SetSecondParticle(FemtoParticle* second){_second = second;}
144+
void SetIdentical(const bool& isidentical){_isidentical = isidentical;}
145+
146+
FemtoParticle* GetFirstParticle() const {return _first;}
147+
FemtoParticle* GetSecondParticle() const {return _second;}
148+
bool IsIdentical(){return _isidentical;}
149+
150+
bool IsClosePair(const float& deta = 0.01, const float& dphi = 0.01, const float& radius = 1.2);
151+
float GetKstar() const;
152+
153+
private:
154+
FemtoParticle* _first;
155+
FemtoParticle* _second;
156+
bool _isidentical = true;
157+
};
158+
159+
bool FemtoPair::IsClosePair(const float& deta, const float& dphi, const float& radius){
160+
if(_first == NULL || _second == NULL) return true;
161+
if(abs(_first->GetEta() - _second->GetEta()) < deta || abs(_first->GetPhiStar(radius) - _second->GetPhiStar(radius)) < dphi) return true;
162+
163+
return false;
164+
}
165+
166+
float FemtoPair::GetKstar() const {
167+
if(_first == NULL || _second == NULL) return -1000;
168+
169+
TLorentzVector* first4momentum = new TLorentzVector( *(_first->Get4momentum()) );
170+
TLorentzVector* second4momentum = new TLorentzVector( *(_second->Get4momentum()) );
171+
172+
if(_isidentical){
173+
TLorentzVector fourmomentadiff = *first4momentum - *second4momentum;
174+
return 0.5*abs(fourmomentadiff.Mag());
175+
}
176+
else{
177+
TLorentzVector fourmomentasum = *first4momentum + *second4momentum;
178+
179+
first4momentum->Boost( (-1)*fourmomentasum.BoostVector() );
180+
second4momentum->Boost( (-1)*fourmomentasum.BoostVector() );
181+
182+
TVector3 qinv = first4momentum->Vect() - second4momentum->Vect();
183+
return 0.5*abs(qinv.Mag());
184+
}
185+
}
186+
187+
#endif // PWGCF_DATAMODEL_SINGLETRACKSELECTOR_H_

‎PWGCF/DataModel/singletrackselector.h‎

Lines changed: 85 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,7 @@
2020
#include "Framework/AnalysisDataModel.h"
2121
#include "Common/DataModel/PIDResponse.h"
2222
#include "Framework/Logger.h"
23+
#include "Common/DataModel/Multiplicity.h"
2324

2425
namespace o2::aod
2526
{
@@ -111,21 +112,26 @@ struct binning {
111112
} // namespace nsigma
112113

113114
DECLARE_SOA_INDEX_COLUMN(Collision, collision); // Index to the collision
114-
DECLARE_SOA_COLUMN(Px, px, float); // Momentum of the track
115-
DECLARE_SOA_COLUMN(Py, py, float); // Momentum of the track
116-
DECLARE_SOA_COLUMN(Pz, pz, float); // Momentum of the track
117-
DECLARE_SOA_COLUMN(P, p, float); // Momentum of the track
118-
DECLARE_SOA_COLUMN(Pt, pt, float); // Momentum of the track
119-
// DECLARE_SOA_COLUMN(PosZ, posZ, float); // vertex position along z
115+
DECLARE_SOA_COLUMN(HasTOF, hasTOF, bool);
116+
DECLARE_SOA_COLUMN(HasITS, hasITS, bool);
117+
DECLARE_SOA_COLUMN(Px, px, float); // Momentum of the track
118+
DECLARE_SOA_COLUMN(Py, py, float); // Momentum of the track
119+
DECLARE_SOA_COLUMN(Pz, pz, float); // Momentum of the track
120+
DECLARE_SOA_COLUMN(P, p, float); // Momentum of the track
121+
DECLARE_SOA_COLUMN(Pt, pt, float); // Momentum of the track
122+
DECLARE_SOA_COLUMN(TPCInnerParam, tpcInnerParam, float); // vertex position along z
123+
DECLARE_SOA_COLUMN(TPCSignal, tpcSignal, float); // vertex position along z
124+
DECLARE_SOA_COLUMN(Beta, beta, float);
120125
DECLARE_SOA_COLUMN(DcaXY, dcaXY, float); // impact parameter of the track
121126
DECLARE_SOA_COLUMN(DcaZ, dcaZ, float); // impact parameter of the track
122127
DECLARE_SOA_COLUMN(TPCNClsFound, tpcNClsFound, float); // Number of TPC clusters
123-
DECLARE_SOA_COLUMN(TPCChi2NCl, tpcChi2NCl, float); // TPC chi2
124-
DECLARE_SOA_COLUMN(ITSNCls, itsNCls, float); // Number of ITS clusters
125-
DECLARE_SOA_COLUMN(ITSChi2NCl, itsChi2NCl, float); // ITS chi2
128+
DECLARE_SOA_COLUMN(TPCCrossedRowsOverFindableCls, tpcCrossedRowsOverFindableCls, float);
129+
DECLARE_SOA_COLUMN(TPCChi2NCl, tpcChi2NCl, float); // TPC chi2
130+
DECLARE_SOA_COLUMN(ITSNCls, itsNCls, float); // Number of ITS clusters
131+
DECLARE_SOA_COLUMN(ITSChi2NCl, itsChi2NCl, float); // ITS chi2
132+
DECLARE_SOA_COLUMN(Sign, sign, short);
126133
DECLARE_SOA_COLUMN(Eta, eta, float);
127134
DECLARE_SOA_COLUMN(Phi, phi, float);
128-
DECLARE_SOA_COLUMN(Sign, sign, float);
129135
DECLARE_SOA_COLUMN(StoredCrossedRows, storedCrossedRows, storedcrossedrows::binning::binned_t);
130136
DECLARE_SOA_COLUMN(StoredTOFNSigmaPr, storedTofNSigmaPr, nsigma::binning::binned_t);
131137
DECLARE_SOA_COLUMN(StoredTPCNSigmaPr, storedTpcNSigmaPr, nsigma::binning::binned_t);
@@ -143,30 +149,80 @@ DECLARE_SOA_DYNAMIC_COLUMN(TOFNSigmaDe, tofNSigmaDe,
143149
DECLARE_SOA_DYNAMIC_COLUMN(TPCNSigmaDe, tpcNSigmaDe,
144150
[](nsigma::binning::binned_t nsigma_binned) -> float { return singletrackselector::unPack<nsigma::binning>(nsigma_binned); });
145151

152+
DECLARE_SOA_DYNAMIC_COLUMN(Energy, energy,
153+
[](float p, float mass) -> float { return sqrt(p * p + mass * mass); });
154+
155+
DECLARE_SOA_DYNAMIC_COLUMN(TrackCuts, trackCuts,
156+
[](float p, float eta, float dcaXY, float dcaZ,
157+
float tpcNClsFound, float tpcChi2NCl, float itsNCls, float tpcCrossedRowsOverFindableCls, storedcrossedrows::binning::binned_t storedCrossedRows,
158+
std::map<std::string, float>* track_cuts) -> bool {
159+
if(p < (*track_cuts)["min_P"] || p > (*track_cuts)["max_P"]) return false;
160+
if(abs(eta) > (*track_cuts)["eta"]) return false;
161+
if(tpcNClsFound < (*track_cuts)["tpcNClsFound"] || tpcChi2NCl > (*track_cuts)["tpcChi2NCl"]) return false;
162+
if(abs(dcaXY) > (*track_cuts)["dcaXY"] || abs(dcaZ) > (*track_cuts)["dcaZ"]) return false;
163+
if(itsNCls < (*track_cuts)["itsNCls"]) return false;
164+
if(singletrackselector::unPackInt<storedcrossedrows::binning>(storedCrossedRows)<(*track_cuts)["crossedrows"]) return false;
165+
if(tpcCrossedRowsOverFindableCls < (*track_cuts)["crossedRows/findableCls"]) return false;
166+
// if(singletrackselector::unPack<storedcrossedrows::binning>(storedCrossedRows)> (*track_cuts)["crossedrows"]) return false;
167+
168+
return true; });
169+
170+
DECLARE_SOA_DYNAMIC_COLUMN(PIDCuts, pidCuts,
171+
[](float pt, float sign, nsigma::binning::binned_t storedTpcNSigmaPr, nsigma::binning::binned_t storedTofNSigmaPr,
172+
nsigma::binning::binned_t storedTpcNSigmaDe, nsigma::binning::binned_t storedTofNSigmaDe,
173+
std::map<std::string, double>* PID_cuts) -> bool {
174+
if(sign != (*PID_cuts)["sign"]) return false;
175+
if(pt < (*PID_cuts)["PIDtrshld"]){
176+
if((*PID_cuts)["particlePDG"] == 2212 && abs(singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaPr)) > (*PID_cuts)["tpcNSigma"]) return false;
177+
if((*PID_cuts)["particlePDG"] == 1000010020 && abs(singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaDe)) > (*PID_cuts)["tpcNSigma"]) return false;
178+
}
179+
else{
180+
if((*PID_cuts)["particlePDG"] == 2212 && sqrt((singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaPr) * singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaPr)) +
181+
(singletrackselector::unPack<nsigma::binning>(storedTofNSigmaPr)*singletrackselector::unPack<nsigma::binning>(storedTofNSigmaPr)))> (*PID_cuts)["tpctofNSigma"]) return false;
182+
//if((*PID_cuts)["particlePDG"] == 2212 && (abs(singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaPr)) > (*PID_cuts)["tpctofNSigma"]
183+
// || abs(singletrackselector::unPack<nsigma::binning>(storedTofNSigmaPr)) > (*PID_cuts)["tpctofNSigma"])) return false;
184+
//if((*PID_cuts)["particlePDG"] == 1000010020 && (abs(singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaDe)) > (*PID_cuts)["tpctofNSigma"]
185+
// || abs(singletrackselector::unPack<nsigma::binning>(storedTofNSigmaDe)) > (*PID_cuts)["tpctofNSigma"])) return false;
186+
if((*PID_cuts)["particlePDG"] == 1000010020 && sqrt((singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaDe) * singletrackselector::unPack<nsigma::binning>(storedTpcNSigmaDe)) +
187+
(singletrackselector::unPack<nsigma::binning>(storedTofNSigmaDe)*singletrackselector::unPack<nsigma::binning>(storedTofNSigmaDe))) > (*PID_cuts)["tpctofNSigma"]) return false;
188+
}
189+
190+
return true; });
191+
146192
// DECLARE_SOA_DYNAMIC_COLUMN(TOFNSigmaPr, tofNSigmaPr,
147193
// [](nsigma::binning::binned_t nsigma_binned) -> float { return nsigma::binning::bin_width * static_cast<float>(nsigma_binned); });
148194
// DECLARE_SOA_DYNAMIC_COLUMN(TPCNSigmaPr, tpcNSigmaPr,
149195
// [](nsigma::binning::binned_t nsigma_binned) -> float { return nsigma::binning::bin_width * static_cast<float>(nsigma_binned); });
150196

197+
DECLARE_SOA_COLUMN(GlobalIndex, globalIndex, int64_t); // Index to the collision
198+
DECLARE_SOA_COLUMN(Mult, mult, int); // Multiplicity of the collision
199+
DECLARE_SOA_COLUMN(PosZ, posZ, int); // Vertex of the collision
200+
151201
} // namespace singletrackselector
152202

153203
DECLARE_SOA_TABLE(SingleTrackSel, "AOD", "STSEL", // Table of the variables for single track selection.
154204
o2::soa::Index<>,
155205
singletrackselector::CollisionId,
206+
singletrackselector::HasITS,
207+
singletrackselector::HasTOF,
156208
singletrackselector::Px,
157209
singletrackselector::Py,
158210
singletrackselector::Pz,
159211
singletrackselector::P,
160212
singletrackselector::Pt,
213+
singletrackselector::TPCInnerParam,
214+
singletrackselector::TPCSignal,
215+
singletrackselector::Beta,
161216
singletrackselector::DcaXY,
162217
singletrackselector::DcaZ,
163218
singletrackselector::TPCNClsFound,
219+
singletrackselector::TPCCrossedRowsOverFindableCls,
164220
singletrackselector::TPCChi2NCl,
165221
singletrackselector::ITSNCls,
166222
singletrackselector::ITSChi2NCl,
223+
singletrackselector::Sign,
167224
singletrackselector::Eta,
168225
singletrackselector::Phi,
169-
singletrackselector::Sign,
170226
singletrackselector::StoredCrossedRows,
171227
singletrackselector::StoredTOFNSigmaPr,
172228
singletrackselector::StoredTPCNSigmaPr,
@@ -176,6 +232,23 @@ DECLARE_SOA_TABLE(SingleTrackSel, "AOD", "STSEL", // Table of the variables for
176232
singletrackselector::TOFNSigmaPr<singletrackselector::StoredTOFNSigmaPr>,
177233
singletrackselector::TPCNSigmaPr<singletrackselector::StoredTPCNSigmaPr>,
178234
singletrackselector::TOFNSigmaDe<singletrackselector::StoredTOFNSigmaDe>,
179-
singletrackselector::TPCNSigmaDe<singletrackselector::StoredTPCNSigmaDe>);
235+
singletrackselector::TPCNSigmaDe<singletrackselector::StoredTPCNSigmaDe>,
236+
singletrackselector::Energy<singletrackselector::P>,
237+
singletrackselector::TrackCuts<singletrackselector::P, singletrackselector::Eta, singletrackselector::DcaXY,
238+
singletrackselector::DcaZ, singletrackselector::TPCNClsFound, singletrackselector::TPCChi2NCl,
239+
singletrackselector::ITSNCls,
240+
singletrackselector::TPCCrossedRowsOverFindableCls,
241+
singletrackselector::StoredCrossedRows>,
242+
singletrackselector::PIDCuts<singletrackselector::Pt, singletrackselector::Sign,
243+
singletrackselector::StoredTPCNSigmaPr,
244+
singletrackselector::StoredTOFNSigmaPr,
245+
singletrackselector::StoredTPCNSigmaDe,
246+
singletrackselector::StoredTOFNSigmaDe>);
247+
248+
DECLARE_SOA_TABLE(SingleCollSel, "AOD", "SCSEL", // Table of the variables for single track selection.
249+
// o2::soa::Index<>,
250+
singletrackselector::GlobalIndex,
251+
singletrackselector::Mult,
252+
singletrackselector::PosZ);
180253
} // namespace o2::aod
181254
#endif // PWGCF_DATAMODEL_SINGLETRACKSELECTOR_H_

0 commit comments

Comments
 (0)