@@ -290,8 +290,11 @@ struct PhotonHBT {
290290 auto photons2_coll = photons2.sliceBy (perCollision2, collision.globalIndex ());
291291 // LOGF(info, "photons1_coll.size() = %d, photons2_coll.size() = %d", photons1_coll.size(), photons2_coll.size());
292292
293- double values_1d[4 ] = {0 .f };
294- double values_3d[8 ] = {0 .f };
293+ float dca_pos1_3d = 999 .f , dca_ele1_3d = 999 .f , dca_ee1_3d = 999 .f ;
294+ float dca_pos2_3d = 999 .f , dca_ele2_3d = 999 .f , dca_ee2_3d = 999 .f ;
295+
296+ double values_1d[6 ] = {0 .f };
297+ double values_3d[10 ] = {0 .f };
295298 if constexpr (pairtype == PairType::kPCMPCM || pairtype == PairType::kPHOSPHOS || pairtype == PairType::kEMCEMC || pairtype == PairType::kDalitzEEDalitzEE ) {
296299 for (auto & cut : cuts1) {
297300 for (auto & paircut : paircuts) {
@@ -300,8 +303,8 @@ struct PhotonHBT {
300303 continue ;
301304 }
302305
303- values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 ;
304- values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 ;
306+ values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 , values_1d[ 2 ] = 0.0 , values_1d[ 3 ] = 0.0 ;
307+ values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 , values_3d[ 2 ] = 0.0 , values_3d[ 3 ] = 0.0 ;
305308 // center-of-mass system (CMS)
306309 ROOT ::Math::PtEtaPhiMVector v1 (g1.pt (), g1.eta (), g1.phi (), 0 .);
307310 ROOT ::Math::PtEtaPhiMVector v2 (g2.pt (), g2.eta (), g2.phi (), 0 .);
@@ -319,15 +322,27 @@ struct PhotonHBT {
319322 values_1d[1 ] = g2.mass ();
320323 values_3d[0 ] = g1.mass ();
321324 values_3d[1 ] = g2.mass ();
325+
326+ dca_pos1_3d = pos1.dca3DinSigma ();
327+ dca_ele1_3d = ele1.dca3DinSigma ();
328+ dca_ee1_3d = std::sqrt ((dca_pos1_3d * dca_pos1_3d + dca_ele1_3d * dca_ele1_3d) / 2 .);
329+ values_1d[2 ] = dca_ee1_3d;
330+ values_3d[2 ] = dca_ee1_3d;
331+
332+ dca_pos2_3d = pos2.dca3DinSigma ();
333+ dca_ele2_3d = ele2.dca3DinSigma ();
334+ dca_ee2_3d = std::sqrt ((dca_pos2_3d * dca_pos2_3d + dca_ele2_3d * dca_ele2_3d) / 2 .);
335+ values_1d[3 ] = dca_ee2_3d;
336+ values_3d[3 ] = dca_ee2_3d;
322337 }
323338
324339 ROOT ::Math::PtEtaPhiMVector q12 = v1 - v2;
325340 ROOT ::Math::PtEtaPhiMVector k12 = 0.5 * (v1 + v2);
326341 float qinv = -q12.M ();
327342 float kt = k12.Pt ();
328343
329- values_1d[2 ] = kt;
330- values_1d[3 ] = qinv;
344+ values_1d[4 ] = kt;
345+ values_1d[5 ] = qinv;
331346
332347 if (fConfigDo3D ) {
333348 // float qt = q12.Pt();
@@ -356,12 +371,12 @@ struct PhotonHBT {
356371 // LOGF(info, "v1.Pz() = %f, v2.Pz() = %f",v1.Pz(), v2.Pz());
357372 // LOGF(info, "v1_lcms_cartesian.Pz() = %f, v2_lcms_cartesian.Pz() = %f",v1_lcms_cartesian.Pz(), v2_lcms_cartesian.Pz());
358373 // LOGF(info, "q12_lcms_cartesian.Pz() = %f", q12_lcms_cartesian.Pz());
359- values_3d[2 ] = kt;
360- values_3d[3 ] = qinv;
361- values_3d[4 ] = qlong_cms;
362- values_3d[5 ] = qout_cms;
363- values_3d[6 ] = qside_cms;
364- values_3d[7 ] = qlong_lcms;
374+ values_3d[4 ] = kt;
375+ values_3d[5 ] = qinv;
376+ values_3d[6 ] = qlong_cms;
377+ values_3d[7 ] = qout_cms;
378+ values_3d[8 ] = qside_cms;
379+ values_3d[9 ] = qlong_lcms;
365380 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut.GetName (), cut.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_same" ))->Fill (values_3d);
366381 } else {
367382 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut.GetName (), cut.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_same" ))->Fill (values_1d);
@@ -388,23 +403,30 @@ struct PhotonHBT {
388403 }
389404 }
390405
391- values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 ;
392- values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 ;
406+ values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 , values_1d[ 2 ] = 0.0 , values_1d[ 3 ] = 0.0 ;
407+ values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 , values_3d[ 2 ] = 0.0 , values_3d[ 3 ] = 0.0 ;
393408 // center-of-mass system (CMS)
394409 ROOT ::Math::PtEtaPhiMVector v1 (g1.pt (), g1.eta (), g1.phi (), 0 .);
395410 ROOT ::Math::PtEtaPhiMVector v2 (g2.pt (), g2.eta (), g2.phi (), 0 .);
396411 if constexpr (pairtype == PairType::kPCMDalitzEE ) {
412+ auto pos2 = g2.template posTrack_as <TEMPrimaryElectrons>();
413+ auto ele2 = g2.template negTrack_as <TEMPrimaryElectrons>();
397414 v2.SetM (g2.mass ());
398415 values_1d[1 ] = g2.mass ();
399416 values_3d[1 ] = g2.mass ();
417+ dca_pos2_3d = pos2.dca3DinSigma ();
418+ dca_ele2_3d = ele2.dca3DinSigma ();
419+ dca_ee2_3d = std::sqrt ((dca_pos2_3d * dca_pos2_3d + dca_ele2_3d * dca_ele2_3d) / 2 .);
420+ values_1d[3 ] = dca_ee2_3d;
421+ values_3d[3 ] = dca_ee2_3d;
400422 }
401423 ROOT ::Math::PtEtaPhiMVector q12 = v1 - v2;
402424 ROOT ::Math::PtEtaPhiMVector k12 = 0.5 * (v1 + v2);
403425 float qinv = -q12.M ();
404426 float kt = k12.Pt ();
405427
406- values_1d[2 ] = kt;
407- values_1d[3 ] = qinv;
428+ values_1d[4 ] = kt;
429+ values_1d[5 ] = qinv;
408430
409431 if (fConfigDo3D ) {
410432 // float qt = q12.Pt();
@@ -425,12 +447,12 @@ struct PhotonHBT {
425447 ROOT ::Math::Boost bst_z (0 , 0 , -beta_z); // Boost supports only PxPyPzEVector
426448 ROOT ::Math::PxPyPzEVector q12_lcms = bst_z (q12_cartesian);
427449 float qlong_lcms = q12_lcms.Pz ();
428- values_3d[2 ] = kt;
429- values_3d[3 ] = qinv;
430- values_3d[4 ] = qlong_cms;
431- values_3d[5 ] = qout_cms;
432- values_3d[6 ] = qside_cms;
433- values_3d[7 ] = qlong_lcms;
450+ values_3d[4 ] = kt;
451+ values_3d[5 ] = qinv;
452+ values_3d[6 ] = qlong_cms;
453+ values_3d[7 ] = qout_cms;
454+ values_3d[8 ] = qside_cms;
455+ values_3d[9 ] = qlong_lcms;
434456 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut1.GetName (), cut2.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_same" ))->Fill (values_3d);
435457 } else {
436458 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut1.GetName (), cut2.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_same" ))->Fill (values_1d);
@@ -478,8 +500,10 @@ struct PhotonHBT {
478500 auto photons_coll2 = photons2.sliceBy (perCollision2, collision2.globalIndex ());
479501 // LOGF(info, "collision1: posZ = %f, numContrib = %d , sel8 = %d | collision2: posZ = %f, numContrib = %d , sel8 = %d", collision1.posZ(), collision1.numContrib(), collision1.sel8(), collision2.posZ(), collision2.numContrib(), collision2.sel8());
480502
481- double values_1d[4 ] = {0 .f };
482- double values_3d[8 ] = {0 .f };
503+ float dca_pos1_3d = 999 .f , dca_ele1_3d = 999 .f , dca_ee1_3d = 999 .f ;
504+ float dca_pos2_3d = 999 .f , dca_ele2_3d = 999 .f , dca_ee2_3d = 999 .f ;
505+ double values_1d[6 ] = {0 .f };
506+ double values_3d[10 ] = {0 .f };
483507 for (auto & cut1 : cuts1) {
484508 for (auto & cut2 : cuts2) {
485509 for (auto & paircut : paircuts) {
@@ -497,28 +521,51 @@ struct PhotonHBT {
497521 }
498522
499523 // center-of-mass system (CMS)
500- values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 ;
501- values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 ;
524+ values_1d[0 ] = 0.0 , values_1d[1 ] = 0.0 , values_1d[ 2 ] = 0.0 , values_1d[ 3 ] = 0.0 ;
525+ values_3d[0 ] = 0.0 , values_3d[1 ] = 0.0 , values_3d[ 2 ] = 0.0 , values_3d[ 3 ] = 0.0 ;
502526 ROOT ::Math::PtEtaPhiMVector v1 (g1.pt (), g1.eta (), g1.phi (), 0 .);
503527 ROOT ::Math::PtEtaPhiMVector v2 (g2.pt (), g2.eta (), g2.phi (), 0 .);
504528 if constexpr (pairtype == PairType::kPCMDalitzEE ) {
529+ auto pos2 = g2.template posTrack_as <TEMPrimaryElectrons>();
530+ auto ele2 = g2.template negTrack_as <TEMPrimaryElectrons>();
505531 v2.SetM (g2.mass ());
506532 values_1d[1 ] = g2.mass ();
507533 values_3d[1 ] = g2.mass ();
534+ dca_pos2_3d = pos2.dca3DinSigma ();
535+ dca_ele2_3d = ele2.dca3DinSigma ();
536+ dca_ee2_3d = std::sqrt ((dca_pos2_3d * dca_pos2_3d + dca_ele2_3d * dca_ele2_3d) / 2 .);
537+ values_1d[3 ] = dca_ee2_3d;
538+ values_3d[3 ] = dca_ee2_3d;
508539 } else if constexpr (pairtype == PairType::kDalitzEEDalitzEE ) {
540+ auto pos1 = g1.template posTrack_as <TEMPrimaryElectrons>();
541+ auto ele1 = g1.template negTrack_as <TEMPrimaryElectrons>();
542+ auto pos2 = g2.template posTrack_as <TEMPrimaryElectrons>();
543+ auto ele2 = g2.template negTrack_as <TEMPrimaryElectrons>();
509544 v1.SetM (g1.mass ());
510545 v2.SetM (g2.mass ());
511546 values_1d[0 ] = g1.mass ();
512547 values_1d[1 ] = g2.mass ();
513548 values_3d[0 ] = g1.mass ();
514549 values_3d[1 ] = g2.mass ();
550+
551+ dca_pos1_3d = pos1.dca3DinSigma ();
552+ dca_ele1_3d = ele1.dca3DinSigma ();
553+ dca_ee1_3d = std::sqrt ((dca_pos1_3d * dca_pos1_3d + dca_ele1_3d * dca_ele1_3d) / 2 .);
554+ values_1d[2 ] = dca_ee1_3d;
555+ values_3d[2 ] = dca_ee1_3d;
556+
557+ dca_pos2_3d = pos2.dca3DinSigma ();
558+ dca_ele2_3d = ele2.dca3DinSigma ();
559+ dca_ee2_3d = std::sqrt ((dca_pos2_3d * dca_pos2_3d + dca_ele2_3d * dca_ele2_3d) / 2 .);
560+ values_1d[3 ] = dca_ee2_3d;
561+ values_3d[3 ] = dca_ee2_3d;
515562 }
516563 ROOT ::Math::PtEtaPhiMVector q12 = v1 - v2;
517564 ROOT ::Math::PtEtaPhiMVector k12 = 0.5 * (v1 + v2);
518565 float qinv = -q12.M ();
519566 float kt = k12.Pt ();
520- values_1d[2 ] = kt;
521- values_1d[3 ] = qinv;
567+ values_1d[4 ] = kt;
568+ values_1d[5 ] = qinv;
522569
523570 if (fConfigDo3D ) {
524571 // float qt = q12.Pt();
@@ -540,12 +587,12 @@ struct PhotonHBT {
540587 ROOT ::Math::PxPyPzEVector q12_lcms = bst_z (q12_cartesian);
541588 float qlong_lcms = q12_lcms.Pz ();
542589
543- values_3d[2 ] = kt;
544- values_3d[3 ] = qinv;
545- values_3d[4 ] = qlong_cms;
546- values_3d[5 ] = qout_cms;
547- values_3d[6 ] = qside_cms;
548- values_3d[7 ] = qlong_lcms;
590+ values_3d[4 ] = kt;
591+ values_3d[5 ] = qinv;
592+ values_3d[6 ] = qlong_cms;
593+ values_3d[7 ] = qout_cms;
594+ values_3d[8 ] = qside_cms;
595+ values_3d[9 ] = qlong_lcms;
549596 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut1.GetName (), cut2.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_mix" ))->Fill (values_3d);
550597 } else {
551598 reinterpret_cast <THnSparseF*>(list_pair_ss->FindObject (Form (" %s_%s" , cut1.GetName (), cut2.GetName ()))->FindObject (paircut.GetName ())->FindObject (" hs_q_mix" ))->Fill (values_1d);
0 commit comments