2424#include " Common/DataModel/EventSelection.h"
2525#include " Common/DataModel/PIDResponse.h"
2626#include " CommonConstants/PhysicsConstants.h"
27+ #include " PWGLF/DataModel/LFKinkDecayTables.h"
2728
2829using namespace o2 ;
2930using namespace o2 ::framework;
3031using namespace o2 ::framework::expressions;
31- using std::array;
32+
33+ using CollisionsFull = soa::Join<aod::Collisions, aod::EvSel>;
3234using FullTracksExtIU = soa::Join<aod::TracksIU, aod::TracksExtra, aod::TracksCovIU, aod::pidTPCFullPr, aod::pidTPCFullAl, aod::pidTPCFullTr, aod::pidTPCFullPi>;
3335using MCLabeledTracksIU = soa::Join<FullTracksExtIU, aod::McTrackLabels>;
3436
@@ -97,17 +99,52 @@ Channel getDecayChannelH4S(TMCParticle const& particle)
9799 return kNDecayChannel ;
98100}
99101// --------------------------------------------------------------
102+ struct Hyperhelium4sigmaAnalysis {
103+ // Histograms are defined with HistogramRegistry
104+ HistogramRegistry registry{" registry" , {}};
105+
106+ // Configurable for event selection
107+ Configurable<float > cutzvertex{" cutzvertex" , 10 .0f , " Accepted z-vertex range (cm)" };
108+ Configurable<float > cutNSigmaAl{" cutNSigmaAl" , 5 , " NSigmaTPCAlpha" };
109+
110+ void init (InitContext const &)
111+ {
112+ // Axes
113+ const AxisSpec vertexZAxis{100 , -15 ., 15 ., " vrtx_{Z} [cm]" };
114+ const AxisSpec ptAxis{50 , -10 , 10 , " #it{p}_{T} (GeV/#it{c})" };
115+ const AxisSpec nSigmaAxis{120 , -6 .f , 6 .f , " n#sigma_{#alpha}" };
116+ const AxisSpec massAxis{100 , 3.85 , 4.25 , " m (GeV/#it{c}^{2})" };
117+
118+ registry.add (" hVertexZRec" , " hVertexZRec" , {HistType::kTH1F , {vertexZAxis}});
100119
120+ registry.add (" h2MassHyperhelium4sigmaPt" , " h2MassHyperhelium4sigmaPt" , {HistType::kTH2F , {ptAxis, massAxis}});
121+ registry.add (" h2NSigmaAlPt" , " h2NSigmaAlPt" , {HistType::kTH2F , {ptAxis, nSigmaAxis}});
122+ }
123+
124+ void process (soa::Join<aod::Collisions, aod::EvSels>::iterator const & collision,
125+ aod::KinkCands const & KinkCands, FullTracksExtIU const &)
126+ {
127+ if (std::abs (collision.posZ ()) > cutzvertex || !collision.sel8 ()) {
128+ return ;
129+ }
130+ registry.fill (HIST (" hVertexZRec" ), collision.posZ ());
131+ for (const auto & kinkCand : KinkCands) {
132+ auto dauTrack = kinkCand.trackDaug_as <FullTracksExtIU>();
133+ if (std::abs (dauTrack.tpcNSigmaAl ()) > cutNSigmaAl) {
134+ continue ;
135+ }
136+ float invMass = RecoDecay::m (std::array{std::array{kinkCand.pxDaug (), kinkCand.pyDaug (), kinkCand.pzDaug ()}, std::array{kinkCand.pxDaugNeut (), kinkCand.pyDaugNeut (), kinkCand.pzDaugNeut ()}}, std::array{o2::constants::physics::MassAlpha, o2::constants::physics::MassPi0});
137+ registry.fill (HIST (" h2MassHyperhelium4sigmaPt" ), kinkCand.mothSign () * kinkCand.ptMoth (), invMass);
138+ registry.fill (HIST (" h2NSigmaAlPt" ), kinkCand.mothSign () * kinkCand.ptDaug (), dauTrack.tpcNSigmaAl ());
139+ }
140+ }
141+ };
142+
143+ // --------------------------------------------------------------
101144// check the performance of mcparticle
102- struct Hyperhelium4sigmaAnalysis {
145+ struct Hyperhelium4sigmaQa {
103146 // Basic checks
104- HistogramRegistry registry{
105- " registry" ,
106- {
107- {" hCollCounter" , " hCollCounter" , {HistType::kTH1F , {{2 , 0 .0f , 2 .0f }}}},
108- {" hMcCollCounter" , " hMcCollCounter" , {HistType::kTH1F , {{2 , 0 .0f , 2 .0f }}}},
109- },
110- };
147+ HistogramRegistry registry{" registry" , {}};
111148
112149 ConfigurableAxis ptBins{" ptBins" , {200 , 0 .f , 10 .f }, " Binning for #it{p}_{T} (GeV/#it{c})" };
113150 ConfigurableAxis ctBins{" ctBins" , {100 , 0 .f , 25 .f }, " Binning for c#it{t} (cm)" };
@@ -117,40 +154,45 @@ struct Hyperhelium4sigmaAnalysis {
117154
118155 void init (InitContext&)
119156 {
120- const AxisSpec ptAxis{ptBins, " #it{p}_{T} (GeV/#it{c})" };
121- const AxisSpec ctAxis{ctBins, " c#it{t} (cm)" };
122- const AxisSpec rigidityAxis{rigidityBins, " p/z (GeV/#it{c})" };
123- const AxisSpec nsigmaAxis{nsigmaBins, " TPC n#sigma" };
124- const AxisSpec invMassAxis{invMassBins, " Inv Mass (GeV/#it{c}^{2})" };
125-
126- registry.get <TH1>(HIST (" hCollCounter" ))->GetXaxis ()->SetBinLabel (1 , " Reconstructed Collisions" );
127- registry.get <TH1>(HIST (" hCollCounter" ))->GetXaxis ()->SetBinLabel (2 , " Selected" );
128- registry.get <TH1>(HIST (" hMcCollCounter" ))->GetXaxis ()->SetBinLabel (1 , " MC Collisions" );
129- registry.get <TH1>(HIST (" hMcCollCounter" ))->GetXaxis ()->SetBinLabel (2 , " Reconstructed" );
130-
131- auto hGenHyperHelium4SigmaCounter = registry.add <TH1>(" hGenHyperHelium4SigmaCounter" , " " , HistType::kTH1F , {{10 , 0 .f , 10 .f }});
132- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (1 , " H4S All" );
133- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (2 , " Matter" );
134- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (3 , " AntiMatter" );
135- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (4 , " #alpha + #pi^{0}" );
136- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (5 , " #bar{#alpha} + #pi^{0}" );
137- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (6 , " t + p + #pi^{0}" );
138- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (7 , " #bar{t} + #bar{p} + #pi^{0}" );
139- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (8 , " t + n + #pi^{+}" );
140- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (9 , " #bar{t} + #bar{n} + #pi^{+}" );
141- registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (10 , " Unexpected" );
142-
143- auto hEvtSelectedHyperHelium4SigmaCounter = registry.add <TH1>(" hEvtSelectedHyperHelium4SigmaCounter" , " " , HistType::kTH1F , {{2 , 0 .f , 2 .f }});
144- registry.get <TH1>(HIST (" hEvtSelectedHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (1 , " Generated" );
145- registry.get <TH1>(HIST (" hEvtSelectedHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (2 , " Survived" );
146-
147- registry.add <TH1>(" hGenHyperHelium4SigmaPt" , " " , HistType::kTH1F , {ptAxis});
148- registry.add <TH1>(" hGenHyperHelium4SigmaCt" , " " , HistType::kTH1F , {ctAxis});
149- registry.add <TH1>(" hMcRecoInvMass" , " " , HistType::kTH1F , {invMassAxis});
150-
151- registry.add <TH2>(" hDauHelium4TPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
152- registry.add <TH2>(" hDauTritonTPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
153- registry.add <TH2>(" hDauProtonTPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
157+ if (doprocessMC == true ) {
158+ const AxisSpec ptAxis{ptBins, " #it{p}_{T} (GeV/#it{c})" };
159+ const AxisSpec ctAxis{ctBins, " c#it{t} (cm)" };
160+ const AxisSpec rigidityAxis{rigidityBins, " p/z (GeV/#it{c})" };
161+ const AxisSpec nsigmaAxis{nsigmaBins, " TPC n#sigma" };
162+ const AxisSpec invMassAxis{invMassBins, " Inv Mass (GeV/#it{c}^{2})" };
163+
164+ auto hCollCounter = registry.add <TH1>(" hCollCounter" , " hCollCounter" , HistType::kTH1F , {{2 , 0 .0f , 2 .0f }});
165+ registry.get <TH1>(HIST (" hCollCounter" ))->GetXaxis ()->SetBinLabel (1 , " Reconstructed Collisions" );
166+ registry.get <TH1>(HIST (" hCollCounter" ))->GetXaxis ()->SetBinLabel (2 , " Selected" );
167+ auto hMcCollCounter = registry.add <TH1>(" hMcCollCounter" , " hMcCollCounter" , HistType::kTH1F , {{2 , 0 .0f , 2 .0f }});
168+ registry.get <TH1>(HIST (" hMcCollCounter" ))->GetXaxis ()->SetBinLabel (1 , " MC Collisions" );
169+ registry.get <TH1>(HIST (" hMcCollCounter" ))->GetXaxis ()->SetBinLabel (2 , " Reconstructed" );
170+
171+ auto hGenHyperHelium4SigmaCounter = registry.add <TH1>(" hGenHyperHelium4SigmaCounter" , " " , HistType::kTH1F , {{10 , 0 .f , 10 .f }});
172+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (1 , " H4S All" );
173+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (2 , " Matter" );
174+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (3 , " AntiMatter" );
175+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (4 , " #alpha + #pi^{0}" );
176+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (5 , " #bar{#alpha} + #pi^{0}" );
177+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (6 , " t + p + #pi^{0}" );
178+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (7 , " #bar{t} + #bar{p} + #pi^{0}" );
179+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (8 , " t + n + #pi^{+}" );
180+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (9 , " #bar{t} + #bar{n} + #pi^{+}" );
181+ registry.get <TH1>(HIST (" hGenHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (10 , " Unexpected" );
182+
183+ auto hEvtSelectedHyperHelium4SigmaCounter = registry.add <TH1>(" hEvtSelectedHyperHelium4SigmaCounter" , " " , HistType::kTH1F , {{2 , 0 .f , 2 .f }});
184+ registry.get <TH1>(HIST (" hEvtSelectedHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (1 , " Generated" );
185+ registry.get <TH1>(HIST (" hEvtSelectedHyperHelium4SigmaCounter" ))->GetXaxis ()->SetBinLabel (2 , " Survived" );
186+
187+ registry.add <TH1>(" hGenHyperHelium4SigmaP" , " " , HistType::kTH1F , {ptAxis});
188+ registry.add <TH1>(" hGenHyperHelium4SigmaPt" , " " , HistType::kTH1F , {ptAxis});
189+ registry.add <TH1>(" hGenHyperHelium4SigmaCt" , " " , HistType::kTH1F , {ctAxis});
190+ registry.add <TH1>(" hMcRecoInvMass" , " " , HistType::kTH1F , {invMassAxis});
191+
192+ registry.add <TH2>(" hDauHelium4TPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
193+ registry.add <TH2>(" hDauTritonTPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
194+ registry.add <TH2>(" hDauProtonTPCNSigma" , " " , HistType::kTH2F , {rigidityAxis, nsigmaAxis});
195+ }
154196 }
155197
156198 // Configurable<bool> eventSel8Cut{"eventSel8Cut", true, "flag to enable event sel8 selection"};
@@ -183,7 +225,13 @@ struct Hyperhelium4sigmaAnalysis {
183225 }
184226 }
185227
186- void process (aod::McCollisions const & mcCollisions, aod::McParticles const & particlesMC, o2::soa::Join<o2::aod::Collisions, o2::aod::McCollisionLabels, o2::aod::EvSels> const & collisions, MCLabeledTracksIU const & tracks)
228+ void processData (o2::aod::Collisions const &)
229+ {
230+ // dummy process function;
231+ }
232+ PROCESS_SWITCH (Hyperhelium4sigmaQa, processData, " process data" , true );
233+
234+ void processMC (aod::McCollisions const & mcCollisions, aod::McParticles const & particlesMC, o2::soa::Join<o2::aod::Collisions, o2::aod::McCollisionLabels, o2::aod::EvSels> const & collisions, MCLabeledTracksIU const & tracks)
187235 {
188236 setTrackIDForMC (particlesMC, tracks);
189237 std::vector<int64_t > selectedEvents (collisions.size ());
@@ -296,6 +344,7 @@ struct Hyperhelium4sigmaAnalysis {
296344 }
297345 }
298346
347+ registry.fill (HIST (" hGenHyperHelium4SigmaP" ), mcparticle.p ());
299348 registry.fill (HIST (" hGenHyperHelium4SigmaPt" ), mcparticle.pt ());
300349 double ct = RecoDecay::sqrtSumOfSquares (svPos[0 ] - mcparticle.vx (), svPos[1 ] - mcparticle.vy (), svPos[2 ] - mcparticle.vz ()) * o2::constants::physics::MassHyperHelium4Sigma / mcparticle.p ();
301350 registry.fill (HIST (" hGenHyperHelium4SigmaCt" ), ct);
@@ -328,11 +377,13 @@ struct Hyperhelium4sigmaAnalysis {
328377 }
329378 }
330379 }
380+ PROCESS_SWITCH (Hyperhelium4sigmaQa, processMC, " do QA for MC prodcutions" , false );
331381};
332382
333383WorkflowSpec defineDataProcessing (ConfigContext const & cfgc)
334384{
335385 return WorkflowSpec{
336386 adaptAnalysisTask<Hyperhelium4sigmaAnalysis>(cfgc),
387+ adaptAnalysisTask<Hyperhelium4sigmaQa>(cfgc),
337388 };
338389}
0 commit comments