diff --git a/PWGCF/JCorran/Tasks/jEPFlowAnalysis.cxx b/PWGCF/JCorran/Tasks/jEPFlowAnalysis.cxx
index b5128bc94eb..be59c91a812 100644
--- a/PWGCF/JCorran/Tasks/jEPFlowAnalysis.cxx
+++ b/PWGCF/JCorran/Tasks/jEPFlowAnalysis.cxx
@@ -65,7 +65,7 @@ struct jEPFlowAnalysis {
 
   Configurable<bool> cfgShiftCorr{"cfgShiftCorr", false, "additional shift correction"};
   Configurable<std::string> cfgShiftPath{"cfgShiftPath", "Users/j/junlee/Qvector/QvecCalib/Shift", "Path for Shift"};
-  Configurable<bool> cfgSPmethod{"cfgSPmethod", false, "flag for scalar product"};
+  Configurable<float> cfgVertexZ{"cfgVertexZ", 10.0, "Maximum vertex Z selection"};
 
   Configurable<std::string> cfgDetName{"cfgDetName", "FT0C", "The name of detector to be analyzed"};
   Configurable<std::string> cfgRefAName{"cfgRefAName", "TPCPos", "The name of detector for reference A"};
@@ -74,6 +74,7 @@ struct jEPFlowAnalysis {
   ConfigurableAxis cfgAxisCent{"cfgAxisCent", {100, 0, 100}, ""};
   ConfigurableAxis cfgAxisPt{"cfgAxisPt", {VARIABLE_WIDTH, 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 10.0, 12.0, 15.0, 30.0, 50.0, 70.0, 100.0}, ""};
   ConfigurableAxis cfgAxisCos{"cfgAxisCos", {102, -1.02, 1.02}, ""};
+  ConfigurableAxis cfgAxisQvec{"cfgAxisQvec", {200, -5.0, 5.0}, ""};
 
   Filter trackFilter = (aod::track::pt > cfgTrackCuts.cfgPtMin) && (nabs(aod::track::eta) < cfgTrackCuts.cfgEtaMax);
 
@@ -123,6 +124,7 @@ struct jEPFlowAnalysis {
     AxisSpec axisCent{cfgAxisCent, "cent"};
     AxisSpec axisPt{cfgAxisPt, "pT"};
     AxisSpec axisCos{cfgAxisCos, "cos"};
+    AxisSpec axisQvec{cfgAxisQvec, "Qvec"};
 
     epFlowHistograms.add("EpDet", "", {HistType::kTH3F, {axisMod, axisCent, axisEvtPl}});
     epFlowHistograms.add("EpRefA", "", {HistType::kTH3F, {axisMod, axisCent, axisEvtPl}});
@@ -135,6 +137,16 @@ struct jEPFlowAnalysis {
     epFlowHistograms.add("vncos", "", {HistType::kTHnSparseF, {axisMod, axisCent, axisPt, axisCos}});
     epFlowHistograms.add("vnsin", "", {HistType::kTHnSparseF, {axisMod, axisCent, axisPt, axisCos}});
 
+    epFlowHistograms.add("EpResQvecDetRefAxx", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+    epFlowHistograms.add("EpResQvecDetRefAxy", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+    epFlowHistograms.add("EpResQvecDetRefBxx", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+    epFlowHistograms.add("EpResQvecDetRefBxy", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+    epFlowHistograms.add("EpResQvecRefARefBxx", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+    epFlowHistograms.add("EpResQvecRefARefBxy", "", {HistType::kTH3F, {axisMod, axisCent, axisQvec}});
+
+    epFlowHistograms.add("SPvnxx", "", {HistType::kTHnSparseF, {axisMod, axisCent, axisPt, axisQvec}});
+    epFlowHistograms.add("SPvnxy", "", {HistType::kTHnSparseF, {axisMod, axisCent, axisPt, axisQvec}});
+
     epFlowHistograms.add("hCentrality", "", {HistType::kTH1F, {axisCent}});
     epFlowHistograms.add("hVertex", "", {HistType::kTH1F, {axisVertex}});
   }
@@ -142,6 +154,8 @@ struct jEPFlowAnalysis {
   void process(MyCollisions::iterator const& coll, soa::Filtered<MyTracks> const& tracks, aod::BCsWithTimestamps const&)
   {
     if (cfgAddEvtSel) {
+      if (std::abs(coll.posZ()) > cfgVertexZ)
+        return;
       switch (cfgEvtSel) {
         case 0: // Sel8
           if (!coll.sel8())
@@ -172,6 +186,8 @@ struct jEPFlowAnalysis {
     epFlowHistograms.fill(HIST("hCentrality"), cent);
     epFlowHistograms.fill(HIST("hVertex"), coll.posZ());
     float eps[3] = {0.};
+    float qx_shifted[3] = {0.};
+    float qy_shifted[3] = {0.};
 
     if (cfgShiftCorr) {
       auto bc = coll.bc_as<aod::BCsWithTimestamps>();
@@ -222,10 +238,14 @@ struct jEPFlowAnalysis {
         eps[0] += deltapsiDet;
         eps[1] += deltapsiRefA;
         eps[2] += deltapsiRefB;
-      }
 
-      if (cfgSPmethod)
-        weight *= std::sqrt(std::pow(coll.qvecRe()[4 * detId + harmInd], 2) + std::pow(coll.qvecIm()[4 * detId + harmInd], 2));
+        qx_shifted[0] = coll.qvecRe()[4 * detId + harmInd] * TMath::Cos(deltapsiDet) - coll.qvecIm()[4 * detId + harmInd] * TMath::Sin(deltapsiDet);
+        qy_shifted[0] = coll.qvecRe()[4 * detId + harmInd] * TMath::Sin(deltapsiDet) + coll.qvecIm()[4 * detId + harmInd] * TMath::Cos(deltapsiDet);
+        qx_shifted[1] = coll.qvecRe()[4 * refAId + harmInd] * TMath::Cos(deltapsiRefA) - coll.qvecIm()[4 * refAId + harmInd] * TMath::Sin(deltapsiRefA);
+        qy_shifted[1] = coll.qvecRe()[4 * refAId + harmInd] * TMath::Sin(deltapsiRefA) + coll.qvecIm()[4 * refAId + harmInd] * TMath::Cos(deltapsiRefA);
+        qx_shifted[2] = coll.qvecRe()[4 * refBId + harmInd] * TMath::Cos(deltapsiRefB) - coll.qvecIm()[4 * refBId + harmInd] * TMath::Sin(deltapsiRefB);
+        qy_shifted[2] = coll.qvecRe()[4 * refBId + harmInd] * TMath::Sin(deltapsiRefB) + coll.qvecIm()[4 * refBId + harmInd] * TMath::Cos(deltapsiRefB);
+      }
 
       float resNumA = helperEP.GetResolution(eps[0], eps[1], i + 2);
       float resNumB = helperEP.GetResolution(eps[0], eps[2], i + 2);
@@ -239,14 +259,22 @@ struct jEPFlowAnalysis {
       epFlowHistograms.fill(HIST("EpResDetRefB"), i + 2, cent, resNumB);
       epFlowHistograms.fill(HIST("EpResRefARefB"), i + 2, cent, resDenom);
 
-      for (int j = 0; j < cfgnMode; j++) { // loop over detectors used
-        for (const auto& track : tracks) {
-          float vn = std::cos((i + 2) * (track.phi() - eps[j]));
-          float vnSin = std::sin((i + 2) * (track.phi() - eps[j]));
+      epFlowHistograms.fill(HIST("EpResQvecDetRefAxx"), i + 2, cent, qx_shifted[0] * qx_shifted[1] + qy_shifted[0] * qy_shifted[1]);
+      epFlowHistograms.fill(HIST("EpResQvecDetRefAxy"), i + 2, cent, qx_shifted[1] * qy_shifted[0] - qx_shifted[0] * qy_shifted[1]);
+      epFlowHistograms.fill(HIST("EpResQvecDetRefBxx"), i + 2, cent, qx_shifted[0] * qx_shifted[2] + qy_shifted[0] * qy_shifted[2]);
+      epFlowHistograms.fill(HIST("EpResQvecDetRefBxy"), i + 2, cent, qx_shifted[2] * qy_shifted[0] - qx_shifted[0] * qy_shifted[2]);
+      epFlowHistograms.fill(HIST("EpResQvecRefARefAxx"), i + 2, cent, qx_shifted[1] * qx_shifted[2] + qy_shifted[1] * qy_shifted[2]);
+      epFlowHistograms.fill(HIST("EpResQvecRefARefAxy"), i + 2, cent, qx_shifted[2] * qy_shifted[1] - qx_shifted[1] * qy_shifted[2]);
 
-          epFlowHistograms.fill(HIST("vncos"), i + 2, cent, track.pt(), vn * weight);
-          epFlowHistograms.fill(HIST("vnsin"), i + 2, cent, track.pt(), vnSin * weight);
-        }
+      for (const auto& track : tracks) {
+        float vn = std::cos((i + 2) * (track.phi() - eps[0]));
+        float vnSin = std::sin((i + 2) * (track.phi() - eps[0]));
+
+        epFlowHistograms.fill(HIST("vncos"), i + 2, cent, track.pt(), vn * weight);
+        epFlowHistograms.fill(HIST("vnsin"), i + 2, cent, track.pt(), vnSin * weight);
+
+        epFlowHistograms.fill(HIST("SPvnxx"), i + 2, cent, track.pt(), (TMath::Cos(track.phi()) * qx_shifted[0] + TMath::Sin(track.phi()) * qy_shifted[0]) * weight);
+        epFlowHistograms.fill(HIST("SPvnxy"), i + 2, cent, track.pt(), (TMath::Sin(track.phi()) * qx_shifted[0] - TMath::Cos(track.phi()) * qy_shifted[0]) * weight);
       }
     }
   }