[Common] Rename ZDCInterCalibTask and update it

Common/DataModel/CMakeLists.txt

@@ -25,6 +25,6 @@ o2physics_add_header_only_library(DataModel
                                           Qvectors.h
                                           MatchMFTFT0.h
                                           MftmchMatchingML.h
-                                          ZDCInterCalib.h
+                                          ZDCExtra.h
                                           EseTable.h
                                           FwdTrackReAlignTables.h)

Common/DataModel/ZDCExtra.h

@@ -0,0 +1,77 @@
+// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
+// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
+// All rights not expressly granted are reserved.
+//
+// This software is distributed under the terms of the GNU General Public
+// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+
+/// \file ZDCExtra.h
+/// \brief ZDC extra table producer
+/// \author Chiara Oppedisano <chiara.oppedisano@cern.ch>, INFN Torino
+/// \author Uliana Dmitrieva <uliana.dmitrieva@cern.ch>, INFN Torino
+
+#ifndef COMMON_DATAMODEL_ZDCEXTRA_H_
+#define COMMON_DATAMODEL_ZDCEXTRA_H_
+
+#include <Framework/AnalysisDataModel.h>
+
+#include <cstdint>
+
+namespace o2::aod
+{
+namespace zdcextra
+{
+DECLARE_SOA_COLUMN(ZnaTowC, znaTowC, float);               //! Common tower ZNA
+DECLARE_SOA_COLUMN(ZnaTow1, znaTow1, float);               //! Tower 1 ZNA
+DECLARE_SOA_COLUMN(ZnaTow2, znaTow2, float);               //! Tower 2 ZNA
+DECLARE_SOA_COLUMN(ZnaTow3, znaTow3, float);               //! Tower 3 ZNA
+DECLARE_SOA_COLUMN(ZnaTow4, znaTow4, float);               //! Tower 4 ZNA
+DECLARE_SOA_COLUMN(ZnaTdc, znaTdc, float);                 //! TDC ZNA
+DECLARE_SOA_COLUMN(ZnaQx, znaQx, float);                   //! Q-vector X ZNA
+DECLARE_SOA_COLUMN(ZnaQy, znaQy, float);                   //! Q-vector Y ZNA
+DECLARE_SOA_COLUMN(ZncTowC, zncTowC, float);               //! Common tower ZNC
+DECLARE_SOA_COLUMN(ZncTow1, zncTow1, float);               //! Tower 1 ZNC
+DECLARE_SOA_COLUMN(ZncTow2, zncTow2, float);               //! Tower 2 ZNC
+DECLARE_SOA_COLUMN(ZncTow3, zncTow3, float);               //! Tower 3 ZNC
+DECLARE_SOA_COLUMN(ZncTow4, zncTow4, float);               //! Tower 4 ZNC
+DECLARE_SOA_COLUMN(ZncTdc, zncTdc, float);                 //! TDC ZNC
+DECLARE_SOA_COLUMN(ZncQx, zncQx, float);                   //! Q-vector X ZNC
+DECLARE_SOA_COLUMN(ZncQy, zncQy, float);                   //! Q-vector Y ZNC
+DECLARE_SOA_COLUMN(Centrality, centrality, float);         //! Centrality
+DECLARE_SOA_COLUMN(Vx, vx, float);                         //! Vertex X
+DECLARE_SOA_COLUMN(Vy, vy, float);                         //! Vertex Y
+DECLARE_SOA_COLUMN(Vz, vz, float);                         //! Vertex Z
+DECLARE_SOA_COLUMN(Timestamp, timestamp, uint64_t);        //! Timestamp
+DECLARE_SOA_COLUMN(SelectionBits, selectionBits, uint8_t); //! Selection Flags
+} // namespace zdcextra
+
+DECLARE_SOA_TABLE(ZdcExtras, "AOD", "ZDCEXTRA", o2::soa::Index<>,
+                  zdcextra::ZnaTowC,
+                  zdcextra::ZnaTow1,
+                  zdcextra::ZnaTow2,
+                  zdcextra::ZnaTow3,
+                  zdcextra::ZnaTow4,
+                  zdcextra::ZnaTdc,
+                  zdcextra::ZnaQx,
+                  zdcextra::ZnaQy,
+                  zdcextra::ZncTowC,
+                  zdcextra::ZncTow1,
+                  zdcextra::ZncTow2,
+                  zdcextra::ZncTow3,
+                  zdcextra::ZncTow4,
+                  zdcextra::ZncTdc,
+                  zdcextra::ZncQx,
+                  zdcextra::ZncQy,
+                  zdcextra::Centrality,
+                  zdcextra::Vx,
+                  zdcextra::Vy,
+                  zdcextra::Vz,
+                  zdcextra::Timestamp,
+                  zdcextra::SelectionBits);
+} // namespace o2::aod
+
+#endif // COMMON_DATAMODEL_ZDCEXTRA_H_

Common/TableProducer/CMakeLists.txt

@@ -45,7 +45,7 @@
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(multcenttable
                     SOURCES multCentTable.cxx
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
@@ -80,12 +80,12 @@
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(propagationservice
                     SOURCES propagationService.cxx
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DCAFitter KFParticle::KFParticle O2Physics::TPCDriftManager
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(propagationservice-run2
                     SOURCES propagationServiceRun2.cxx
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DCAFitter KFParticle::KFParticle O2Physics::TPCDriftManager
                     COMPONENT_NAME Analysis)
@@ -95,7 +95,7 @@
                     PUBLIC_LINK_LIBRARIES O2::DetectorsBase O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(calo-clusters
                     SOURCES caloClusterProducer.cxx
                     PUBLIC_LINK_LIBRARIES O2::DataFormatsPHOS O2::PHOSBase O2::PHOSReconstruction O2Physics::DataModel
                     COMPONENT_NAME Analysis)
@@ -118,12 +118,12 @@
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(mccollisionextra
                     SOURCES mcCollsExtra.cxx
                     PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(qvector-table
                     SOURCES qVectorsTable.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework
                                           O2Physics::AnalysisCore
@@ -132,7 +132,7 @@
                                           O2::CCDB
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(mftmchmatchingml
                     SOURCES mftmchMatchingML.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework
                                           O2Physics::AnalysisCore
@@ -141,15 +141,15 @@
                                           O2::ReconstructionDataFormats
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(match-mft-ft0
                     SOURCES match-mft-ft0.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
                                           O2::ReconstructionDataFormats
                                           O2::DetectorsBase O2::DetectorsCommonDataFormats
                     COMPONENT_NAME Analysis)
 
-o2physics_add_dpl_workflow(zdc-task-inter-calib
-                    SOURCES zdcTaskInterCalib.cxx
+o2physics_add_dpl_workflow(zdc-extra-table-producer
+                    SOURCES zdcExtraTableProducer.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
@@ -158,12 +158,12 @@
                     PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(mftmch-matching-data
                     SOURCES match-mft-mch-data.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2::GlobalTracking
                     COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(occ-table-producer
                     SOURCES occupancyTableProducer.cxx
                     PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
                     COMPONENT_NAME Analysis)

Common/TableProducer/zdcTaskInterCalib.cxx → Common/TableProducer/zdcExtraTableProducer.cxx

@@ -9,14 +9,15 @@
 // granted to it by virtue of its status as an Intergovernmental Organization
 // or submit itself to any jurisdiction.
 
-/// \file zdcTaskInterCalib.cxx
-/// \brief Task for ZDC tower inter-calibration
-/// \author chiara.oppedisano@cern.ch
+/// \file zdcExtraTableProducer.cxx
+/// \brief Task creating table with ZDC PMTs energies and calculated centroid (Q-vector) to be used for spectator plane measurement
+/// \author Chiara Oppedisano <chiara.oppedisano@cern.ch>, INFN Torino
+/// \author Uliana Dmitrieva <uliana.dmitrieva@cern.ch>, INFN Torino
 
 #include "Common/CCDB/EventSelectionParams.h"
 #include "Common/DataModel/Centrality.h"
 #include "Common/DataModel/EventSelection.h"
-#include "Common/DataModel/ZDCInterCalib.h"
+#include "Common/DataModel/ZDCExtra.h"
 
 #include <Framework/AnalysisDataModel.h>
 #include <Framework/AnalysisHelpers.h>
@@ -29,6 +30,7 @@
 #include <Framework/runDataProcessing.h>
 
 #include <TH1.h>
+#include <TH2.h>
 
 #include <cstdint>
 
@@ -40,17 +42,18 @@ using namespace o2::aod::evsel;
 using BCsRun3 = soa::Join<aod::BCs, aod::Timestamps, aod::BcSels, aod::Run3MatchedToBCSparse>;
 using ColEvSels = soa::Join<aod::Collisions, aod::EvSels, aod::CentFT0Cs>;
 
-struct ZdcTaskInterCalib {
+struct ZdcExtraTableProducer {
 
-  Produces<aod::ZDCInterCalib> zTab;
+  Produces<aod::ZdcExtras> zdcextras;
 
   // Configurable parameters
   //
   Configurable<int> nBins{"nBins", 400, "n bins"};
   Configurable<float> maxZN{"maxZN", 399.5, "Max ZN signal"};
   Configurable<bool> tdcCut{"tdcCut", false, "Flag for TDC cut"};
   Configurable<float> tdcZNmincut{"tdcZNmincut", -2.5, "Min ZN TDC cut"};
-  Configurable<float> tdcZNmaxcut{"tdcZNmaxcut", -2.5, "Max ZN TDC cut"};
+  Configurable<float> tdcZNmaxcut{"tdcZNmaxcut", 2.5, "Max ZN TDC cut"};
+  Configurable<bool> cfgUsePMC{"cfgUsePMC", true, "Use common PM (true) or sum of PMs (false) "};
   // Event selections
   Configurable<bool> cfgEvSelSel8{"cfgEvSelSel8", true, "Event selection: sel8"};
   Configurable<float> cfgEvSelVtxZ{"cfgEvSelVtxZ", 10, "Event selection: zVtx"};
@@ -92,6 +95,9 @@ struct ZdcTaskInterCalib {
     registry.add("ZNAsumq", "ZNAsumq; ZNA uncalib. sum PMQ; Entries", {HistType::kTH1F, {{nBins, -0.5, maxZN}}});
     registry.add("ZNCsumq", "ZNCsumq; ZNC uncalib. sum PMQ; Entries", {HistType::kTH1F, {{nBins, -0.5, maxZN}}});
 
+    registry.add("ZNACentroid", "ZNACentroid; ZNA Centroid; X; Y", {HistType::kTH2F, {{50, -1.5, 1.5}, {50, -1.5, 1.5}}});
+    registry.add("ZNCCentroid", "ZNCCentroid; ZNC Centroid; X; Y", {HistType::kTH2F, {{50, -1.5, 1.5}, {50, -1.5, 1.5}}});
+
     registry.add("hEventCount", "Number of Event; Cut; #Events Passed Cut", {HistType::kTH1D, {{nEventSelections, 0, nEventSelections}}});
     registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_allEvents + 1, "All events");
     registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_zvtx + 1, "vtxZ");
@@ -101,7 +107,7 @@ struct ZdcTaskInterCalib {
     registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_kIsGoodZvtxFT0vsPV + 1, "kIsGoodZvtxFT0vsPV");
     registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_kNoCollInTimeRangeStandard + 1, "kNoCollInTimeRangeStandard");
     registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_kIsVertexITSTPC + 1, "kIsVertexITSTPC");
-    registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_kIsGoodITSLayersAll + 1, "kkIsGoodITSLayersAll");
+    registry.get<TH1>(HIST("hEventCount"))->GetXaxis()->SetBinLabel(evSel_kIsGoodITSLayersAll + 1, "kIsGoodITSLayersAll");
   }
 
   template <typename TCollision>
@@ -206,18 +212,8 @@ struct ZdcTaskInterCalib {
         //
         double sumZNC = 0;
         double sumZNA = 0;
-        double pmqZNC[4] = {
-          0,
-          0,
-          0,
-          0,
-        };
-        double pmqZNA[4] = {
-          0,
-          0,
-          0,
-          0,
-        };
+        double pmqZNC[4] = {};
+        double pmqZNA[4] = {};
         //
         if (isZNChit) {
           for (int it = 0; it < nTowers; it++) {
@@ -244,15 +240,87 @@ struct ZdcTaskInterCalib {
           registry.get<TH1>(HIST("ZNApm4"))->Fill(pmqZNA[3]);
           registry.get<TH1>(HIST("ZNAsumq"))->Fill(sumZNA);
         }
-        if (isZNAhit || isZNChit)
-          zTab(pmcZNA, pmqZNA[0], pmqZNA[1], pmqZNA[2], pmqZNA[3], tdcZNC, pmcZNC, pmqZNC[0], pmqZNC[1], pmqZNC[2], pmqZNC[3], tdcZNA, centrality, foundBC.timestamp(), evSelection);
+
+        // Q-vectors (centroid) calculation
+        // kBeamEne --  LHC Run 3 Pb-Pb collision energy (5.36 TeV per nucleon pair)
+        constexpr float kBeamEne = 5.36 * 0.5;
+
+        // Provide coordinates of centroid over ZN (side C) front face
+        constexpr float X[4] = {-1.75, 1.75, -1.75, 1.75};
+        constexpr float Y[4] = {-1.75, -1.75, 1.75, 1.75};
+        constexpr float kAlpha = 0.395; // saturation correction
+
+        float numXZNC = 0., numYZNC = 0., denZNC = 0.;
+        float numXZNA = 0., numYZNA = 0., denZNA = 0.;
+
+        // Calculate weighted sums of the x and y coordinates
+        constexpr int kNTowers = 4; // number of ZDC towers
+        for (int i = 0; i < kNTowers; i++) {
+          if (pmqZNC[i] > 0.) {
+            float wZNC = std::pow(pmqZNC[i], kAlpha);
+            numXZNC -= X[i] * wZNC; // numerator x (minus sign due to opposite orientation of ZNC)
+            numYZNC += Y[i] * wZNC; // numerator y
+            denZNC += wZNC;         // denominator
+          }
+          if (pmqZNA[i] > 0.) {
+            float wZNA = std::pow(pmqZNA[i], kAlpha);
+            numXZNA += X[i] * wZNA; // numerator x
+            numYZNA += Y[i] * wZNA; // numerator y
+            denZNA += wZNA;         // denominator
+          }
+        }
+        // Calculate centroid coordinates (in cm) with correction factor c depending on the number of spectator nucleons (nSpec)
+
+        float zncCommon = 0;
+        float znaCommon = 0;
+
+        // Use sum of PMTs (cfgUsePMC == false) when common PMT is saturated
+        if (cfgUsePMC) {
+          zncCommon = pmcZNC;
+          znaCommon = pmcZNA;
+        } else {
+          zncCommon = sumZNC;
+          znaCommon = sumZNA;
+        }
+
+        float centroidZNC[2], centroidZNA[2];
+
+        if (denZNC != 0.) {
+          float nSpecnC = zncCommon / kBeamEne;
+          float cZNC = 1.89358 - 0.71262 / (nSpecnC + 0.71789);
+          centroidZNC[0] = cZNC * numXZNC / denZNC;
+          centroidZNC[1] = cZNC * numYZNC / denZNC;
+        } else {
+          centroidZNC[0] = 999.;
+          centroidZNC[1] = 999.;
+        }
+        //
+        if (denZNA != 0.) {
+          float nSpecnA = znaCommon / kBeamEne;
+          float cZNA = 1.89358 - 0.71262 / (nSpecnA + 0.71789);
+          centroidZNA[0] = cZNA * numXZNA / denZNA;
+          centroidZNA[1] = cZNA * numYZNA / denZNA;
+        } else {
+          centroidZNA[0] = 999.;
+          centroidZNA[1] = 999.;
+        }
+        registry.get<TH2>(HIST("ZNCCentroid"))->Fill(centroidZNC[0], centroidZNC[1]);
+        registry.get<TH2>(HIST("ZNACentroid"))->Fill(centroidZNA[0], centroidZNA[1]);
+
+        auto vz = collision.posZ();
+        auto vx = collision.posX();
+        auto vy = collision.posY();
+
+        if (isZNAhit || isZNChit) {
+          zdcextras(pmcZNA, pmqZNA[0], pmqZNA[1], pmqZNA[2], pmqZNA[3], tdcZNA, centroidZNA[0], centroidZNA[1], pmcZNC, pmqZNC[0], pmqZNC[1], pmqZNC[2], pmqZNC[3], tdcZNC, centroidZNC[0], centroidZNC[1], centrality, vx, vy, vz, foundBC.timestamp(), evSelection);
+        }
       }
     }
   }
 };
 
-WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) // o2-linter: disable=name/file-cpp
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
 {
   return WorkflowSpec{
-    adaptAnalysisTask<ZdcTaskInterCalib>(cfgc)};
+    adaptAnalysisTask<ZdcExtraTableProducer>(cfgc)};
 }