99// granted to it by virtue of its status as an Intergovernmental Organization
1010// or submit itself to any jurisdiction.
1111
12- // author: Noor Koster [email protected] 12+ // / \file flowSP.cxx
13+ // / \author Noor Koster
14+ // / \since 01/12/2024
15+ // / \brief task to evaluate flow with respect to spectator plane.
1316
1417#include < CCDB/BasicCCDBManager.h>
1518#include < DataFormatsParameters/GRPObject.h>
3033#include " Common/DataModel/Multiplicity.h"
3134#include " Common/DataModel/Centrality.h"
3235#include " Common/DataModel/Qvectors.h"
36+ #include " Common/Core/RecoDecay.h"
3337
3438#include " PWGCF/DataModel/SPTableZDC.h"
3539#include " TF1.h"
@@ -41,7 +45,7 @@ using namespace o2::framework::expressions;
4145
4246#define O2_DEFINE_CONFIGURABLE (NAME, TYPE, DEFAULT, HELP ) Configurable<TYPE> NAME{#NAME, DEFAULT, HELP};
4347
44- struct flowAnalysisSP {
48+ struct FlowSP {
4549
4650 O2_DEFINE_CONFIGURABLE (cfgDCAxy, float , 0.2 , " Cut on DCA in the transverse direction (cm)" );
4751 O2_DEFINE_CONFIGURABLE (cfgDCAz, float , 2 , " Cut on DCA in the longitudinal direction (cm)" );
@@ -53,6 +57,8 @@ struct flowAnalysisSP {
5357 O2_DEFINE_CONFIGURABLE (cfgMagField, float , 99999 , " Configurable magnetic field; default CCDB will be queried" );
5458 O2_DEFINE_CONFIGURABLE (cfgUseAdditionalEventCut, bool , true , " Bool to enable Additional Event Cut" );
5559 O2_DEFINE_CONFIGURABLE (cfgUseAdditionalTrackCut, bool , true , " Bool to enable Additional Track Cut" );
60+ O2_DEFINE_CONFIGURABLE (cfgCentMax, float , 60 , " Maximum cenrality for selected events" );
61+ O2_DEFINE_CONFIGURABLE (cfgCentMin, float , 10 , " Minimum cenrality for selected events" );
5662
5763 O2_DEFINE_CONFIGURABLE (cfgDoubleTrackFunction, bool , true , " Include track cut at low pt" );
5864 O2_DEFINE_CONFIGURABLE (cfgTrackCutSize, float , 0.06 , " Spread of track cut" );
@@ -67,8 +73,8 @@ struct flowAnalysisSP {
6773
6874 Filter collisionFilter = nabs(aod::collision::posZ) < cfgVtxZ;
6975 Filter trackFilter = nabs(aod::track::eta) < cfgEta && aod::track::pt > cfgPtmin&& aod::track::pt < cfgPtmax && ((requireGlobalTrackInFilter()) || (aod::track::isGlobalTrackSDD == (uint8_t ) true )) && nabs(aod::track::dcaXY) < cfgDCAxy&& nabs(aod::track::dcaZ) < cfgDCAz;
70- using myCollisions = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::Mults, aod::CentFT0Cs, aod::SPTableZDC, aod::Qvectors>>;
71- using myTracks = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TrackSelection, aod::TracksDCA>>;
76+ using UsedCollisions = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::Mults, aod::CentFT0Cs, aod::SPTableZDC, aod::Qvectors>>;
77+ using UsedTracks = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TrackSelection, aod::TracksDCA>>;
7278
7379 // Connect to ccdb
7480 Service<ccdb::BasicCCDBManager> ccdb;
@@ -118,6 +124,22 @@ struct flowAnalysisSP {
118124 registry.add <TProfile>(" v1C_eta" , " " , kTProfile , {{10 , -.8 , .8 }});
119125 registry.add <TProfile>(" v1AC_eta" , " " , kTProfile , {{10 , -.8 , .8 }});
120126
127+ registry.add <TProfile>(" v1_eta_odd" , " " , kTProfile , {{10 , -.8 , .8 }});
128+ registry.add <TProfile>(" v1_eta_even" , " " , kTProfile , {{10 , -.8 , .8 }});
129+
130+ registry.add <TProfile>(" v1_eta_odd_dev" , " " , kTProfile , {{10 , -.8 , .8 }});
131+ registry.add <TProfile>(" v1_eta_even_dev" , " " , kTProfile , {{10 , -.8 , .8 }});
132+
133+ registry.add <TProfile>(" v1_eta_odd_dev_pos" , " " , kTProfile , {{10 , -.8 , .8 }});
134+ registry.add <TProfile>(" v1_eta_even_dev_pos" , " " , kTProfile , {{10 , -.8 , .8 }});
135+ registry.add <TProfile>(" v1_eta_odd_pos" , " " , kTProfile , {{10 , -.8 , .8 }});
136+ registry.add <TProfile>(" v1_eta_even_pos" , " " , kTProfile , {{10 , -.8 , .8 }});
137+
138+ registry.add <TProfile>(" v1_eta_odd_dev_neg" , " " , kTProfile , {{10 , -.8 , .8 }});
139+ registry.add <TProfile>(" v1_eta_even_dev_neg" , " " , kTProfile , {{10 , -.8 , .8 }});
140+ registry.add <TProfile>(" v1_eta_odd_neg" , " " , kTProfile , {{10 , -.8 , .8 }});
141+ registry.add <TProfile>(" v1_eta_even_neg" , " " , kTProfile , {{10 , -.8 , .8 }});
142+
121143 registry.add <TProfile>(" v2_cent" , " " , kTProfile , {{80 , 0 , 80 }});
122144 registry.add <TProfile>(" v2A_cent" , " " , kTProfile , {{80 , 0 , 80 }});
123145 registry.add <TProfile>(" v2C_cent" , " " , kTProfile , {{80 , 0 , 80 }});
@@ -133,7 +155,8 @@ struct flowAnalysisSP {
133155 registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (7 , " kNoCollInTimeRangeStandard" );
134156 registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (8 , " kIsVertexITSTPC" );
135157 registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (9 , " after Mult cuts" );
136- registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (10 , " isSelected" );
158+ registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (10 , " after Cent cuts" );
159+ registry.get <TH1>(HIST (" hEventCount" ))->GetXaxis ()->SetBinLabel (11 , " isSelected" );
137160
138161 if (cfgUseAdditionalEventCut) {
139162 fMultPVCutLow = new TF1 (" fMultPVCutLow" , " [0]+[1]*x+[2]*x*x+[3]*x*x*x+[4]*x*x*x*x - 3.5*([5]+[6]*x+[7]*x*x+[8]*x*x*x+[9]*x*x*x*x)" , 0 , 100 );
@@ -217,7 +240,7 @@ struct flowAnalysisSP {
217240 float vtxz = -999 ;
218241 if (collision.numContrib () > 1 ) {
219242 vtxz = collision.posZ ();
220- float zRes = TMath::Sqrt (collision.covZZ ());
243+ float zRes = std::sqrt (collision.covZZ ());
221244 if (zRes > 0.25 && collision.numContrib () < 20 )
222245 vtxz = -999 ;
223246 }
@@ -239,6 +262,11 @@ struct flowAnalysisSP {
239262
240263 registry.fill (HIST (" hEventCount" ), 8.5 );
241264
265+ if (centrality > cfgCentMax || centrality < cfgCentMin)
266+ return 0 ;
267+
268+ registry.fill (HIST (" hEventCount" ), 9.5 );
269+
242270 return 1 ;
243271 }
244272
@@ -247,22 +275,22 @@ struct flowAnalysisSP {
247275 {
248276 double phimodn = track.phi ();
249277 if (field < 0 ) // for negative polarity field
250- phimodn = TMath::TwoPi () - phimodn;
278+ phimodn = o2::constants::math::TwoPI - phimodn;
251279 if (track.sign () < 0 ) // for negative charge
252- phimodn = TMath::TwoPi () - phimodn;
280+ phimodn = o2::constants::math::TwoPI - phimodn;
253281 if (phimodn < 0 )
254282 LOGF (warning, " phi < 0: %g" , phimodn);
255283
256- phimodn += TMath::Pi () / 18.0 ; // to center gap in the middle
257- phimodn = fmod (phimodn, TMath::Pi () / 9.0 );
284+ phimodn += o2::constants::math::PI / 18.0 ; // to center gap in the middle
285+ phimodn = fmod (phimodn, o2::constants::math::PI / 9.0 );
258286 registry.fill (HIST (" pt_phi_bef" ), track.pt (), phimodn);
259287 if (phimodn < fPhiCutHigh ->Eval (track.pt ()) && phimodn > fPhiCutLow ->Eval (track.pt ()))
260288 return false ; // reject track
261289 registry.fill (HIST (" pt_phi_aft" ), track.pt (), phimodn);
262290 return true ;
263291 }
264292
265- void process (myCollisions ::iterator const & collision, aod::BCsWithTimestamps const &, myTracks const & tracks)
293+ void process (UsedCollisions ::iterator const & collision, aod::BCsWithTimestamps const &, UsedTracks const & tracks)
266294 {
267295 // Hier sum over collisions and get ZDC data.
268296 registry.fill (HIST (" hEventCount" ), .5 );
@@ -280,7 +308,7 @@ struct flowAnalysisSP {
280308 return ;
281309
282310 if (collision.isSelected ()) {
283- registry.fill (HIST (" hEventCount" ), 9 .5 );
311+ registry.fill (HIST (" hEventCount" ), 10 .5 );
284312
285313 registry.fill (HIST (" hCent" ), centrality);
286314
@@ -289,38 +317,90 @@ struct flowAnalysisSP {
289317 double qxC = collision.qxC ();
290318 double qyC = collision.qyC ();
291319
292- double Psi_A = 1.0 * TMath::ATan2 (qyA, qxA);
293- registry.fill (HIST (" hSPplaneA" ), Psi_A, 1 );
320+ double psiA = 1.0 * std::atan2 (qyA, qxA);
321+ registry.fill (HIST (" hSPplaneA" ), psiA, 1 );
322+
323+ double psiC = 1.0 * std::atan2 (qyC, qxC);
324+ registry.fill (HIST (" hSPplaneC" ), psiC, 1 );
325+
326+ registry.fill (HIST (" hSPplaneA-C" ), psiA - psiC, 1 );
294327
295- double Psi_C = 1.0 * TMath::ATan2 (qyC, qxC);
296- registry.fill (HIST (" hSPplaneC" ), Psi_C, 1 );
328+ registry.fill (HIST (" hCosdPhi" ), centrality, std::cos (psiA - psiC));
329+ if (std::cos (psiA - psiC) < 0 )
330+ registry.fill (HIST (" hSPlaneRes" ), centrality, std::sqrt (-1 . * std::cos (psiA - psiC)));
297331
298- registry.fill (HIST (" hSPplaneA-C " ), Psi_A - Psi_C, 1 );
332+ registry.fill (HIST (" hSindPhi " ), centrality, std::sin (psiA - psiC) );
299333
300- registry.fill (HIST (" hCosdPhi" ), centrality, TMath::Cos (Psi_A - Psi_C));
301- if (TMath::Cos (Psi_A - Psi_C) < 0 )
302- registry.fill (HIST (" hSPlaneRes" ), centrality, TMath::Sqrt (-1 . * TMath::Cos (Psi_A - Psi_C)));
303- registry.fill (HIST (" hSindPhi" ), centrality, TMath::Sin (Psi_A - Psi_C));
334+ auto qxAqxC = qxA * qxC;
335+ auto qyAqyC = qyA * qyC;
304336
305- for (auto & track : tracks) {
337+ for (const auto & track : tracks) {
306338 if (!trackSelected (track, field))
307339 continue ;
308340
309- double v1A = TMath::Cos (track.phi () - Psi_A);
310- double v1C = TMath::Cos (track.phi () - Psi_C);
341+ bool pos;
342+ if (track.sign () == 0.0 )
343+ continue ;
344+ if (track.sign () > 0 ) {
345+ pos = true ;
346+ } else {
347+ pos = false ;
348+ }
349+
350+ // constrain angle to 0 -> [0,0+2pi]
351+ auto phi = RecoDecay::constrainAngle (track.phi (), 0 );
311352
312- double v1AC = TMath::Cos (track.phi () - (Psi_A - Psi_C));
353+ auto ux = std::cos (phi);
354+ auto uy = std::sin (phi);
313355
314- registry.fill (HIST (" v1_eta" ), track.eta (), (1 . / TMath::Sqrt (2 )) * (v1A - v1C));
356+ auto uxQxA = ux * qxA;
357+ auto uyQyA = uy * qyA;
358+ auto uxyQxyA = uxQxA + uyQyA;
359+ auto uxQxC = ux * qxC;
360+ auto uyQyC = uy * qyC;
361+ auto uxyQxyC = uxQxC + uyQyC;
362+
363+ auto oddv1 = ux * (qxA - qxC) + uy * (qyA - qyC);
364+ auto evenv1 = ux * (qxA + qxC) + uy * (qyA + qyC);
365+
366+ auto oddv1Dev = ux * (qxA - qxC) / std::sqrt (std::abs (qxAqxC)) + uy * (qyA - qyC) / std::sqrt (std::abs (qyAqyC));
367+ auto evenv1Dev = ux * (qxA + qxC) / std::sqrt (std::abs (qxAqxC)) + uy * (qyA + qyC) / std::sqrt (std::abs (qyAqyC));
368+
369+ double v1A = std::cos (phi - psiA);
370+ double v1C = std::cos (phi - psiC);
371+
372+ double v1AC = std::cos (phi - (psiA - psiC));
373+
374+ registry.fill (HIST (" v1_eta" ), track.eta (), (1 . / std::sqrt (2 )) * (v1A - v1C));
315375 registry.fill (HIST (" v1A_eta" ), track.eta (), (v1A));
316376 registry.fill (HIST (" v1C_eta" ), track.eta (), (v1C));
317377 registry.fill (HIST (" v1AC_eta" ), track.eta (), (v1AC));
318378
319- double v2A = TMath::Cos (2 * (track.phi () - Psi_A));
320- double v2C = TMath::Cos (2 * (track.phi () - Psi_C));
321- double v2AC = TMath::Cos (2 * (track.phi () - (Psi_A - Psi_C)));
379+ registry.fill (HIST (" v1_eta_odd" ), track.eta (), oddv1);
380+ registry.fill (HIST (" v1_eta_even" ), track.eta (), evenv1);
381+
382+ registry.fill (HIST (" v1_eta_odd_dev" ), track.eta (), oddv1Dev);
383+ registry.fill (HIST (" v1_eta_even_dev" ), track.eta (), evenv1Dev);
384+
385+ if (pos) {
386+ registry.fill (HIST (" v1_eta_odd_pos" ), track.eta (), oddv1);
387+ registry.fill (HIST (" v1_eta_even_pos" ), track.eta (), evenv1);
388+
389+ registry.fill (HIST (" v1_eta_odd_dev_pos" ), track.eta (), oddv1Dev);
390+ registry.fill (HIST (" v1_eta_even_dev_pos" ), track.eta (), evenv1Dev);
391+ } else {
392+ registry.fill (HIST (" v1_eta_odd_neg" ), track.eta (), oddv1);
393+ registry.fill (HIST (" v1_eta_even_neg" ), track.eta (), evenv1);
394+
395+ registry.fill (HIST (" v1_eta_odd_dev_neg" ), track.eta (), oddv1Dev);
396+ registry.fill (HIST (" v1_eta_even_dev_neg" ), track.eta (), evenv1Dev);
397+ }
398+
399+ double v2A = std::cos (2 * (phi - psiA));
400+ double v2C = std::cos (2 * (phi - psiC));
401+ double v2AC = std::cos (2 * (phi - (psiA - psiC)));
322402
323- registry.fill (HIST (" v2_cent" ), centrality, (1 . / TMath::Sqrt (2 )) * (v2A - v2C));
403+ registry.fill (HIST (" v2_cent" ), centrality, (1 . / std::sqrt (2 )) * (v2A - v2C));
324404 registry.fill (HIST (" v2A_cent" ), centrality, (v2A));
325405 registry.fill (HIST (" v2C_cent" ), centrality, (v2C));
326406 registry.fill (HIST (" v2AC_cent" ), centrality, (v2AC));
@@ -338,6 +418,6 @@ struct flowAnalysisSP {
338418WorkflowSpec defineDataProcessing (ConfigContext const & cfgc)
339419{
340420 return WorkflowSpec{
341- adaptAnalysisTask<flowAnalysisSP >(cfgc),
421+ adaptAnalysisTask<FlowSP >(cfgc),
342422 };
343423}
0 commit comments