|
| 1 | +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. |
| 2 | +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. |
| 3 | +// All rights not expressly granted are reserved. |
| 4 | +// |
| 5 | +// This software is distributed under the terms of the GNU General Public |
| 6 | +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". |
| 7 | +// |
| 8 | +// In applying this license CERN does not waive the privileges and immunities |
| 9 | +// granted to it by virtue of its status as an Intergovernmental Organization |
| 10 | +// or submit itself to any jurisdiction. |
| 11 | + |
| 12 | +// jet v2 task |
| 13 | +/// \author Yubiao Wang <[email protected]> |
| 14 | +// C++/ROOT includes. |
| 15 | + |
| 16 | +#include <chrono> |
| 17 | +#include <string> |
| 18 | +#include <vector> |
| 19 | +#include <TComplex.h> |
| 20 | +#include <TH3F.h> |
| 21 | +#include <cmath> |
| 22 | +#include <TRandom3.h> |
| 23 | +#include <TF1.h> |
| 24 | +#include <algorithm> |
| 25 | +#include <unordered_map> |
| 26 | +// o2Physics includes. |
| 27 | + |
| 28 | +#include "CCDB/BasicCCDBManager.h" |
| 29 | +#include "DataFormatsParameters/GRPMagField.h" |
| 30 | + |
| 31 | +#include "Framework/runDataProcessing.h" |
| 32 | + |
| 33 | +#include "Common/DataModel/FT0Corrected.h" |
| 34 | +#include "Common/DataModel/Multiplicity.h" |
| 35 | +#include "Common/DataModel/Centrality.h" |
| 36 | +#include "Common/CCDB/ctpRateFetcher.h" |
| 37 | + |
| 38 | +//< evt pln .h >// |
| 39 | +#include "Framework/ASoAHelpers.h" |
| 40 | +#include "Framework/RunningWorkflowInfo.h" |
| 41 | +#include "Framework/StaticFor.h" |
| 42 | + |
| 43 | +#include "Common/DataModel/Qvectors.h" |
| 44 | +#include "Common/Core/EventPlaneHelper.h" |
| 45 | +//< evt pln .h | end >// |
| 46 | + |
| 47 | +// o2 includes. |
| 48 | +#include "DetectorsCommonDataFormats/AlignParam.h" |
| 49 | + |
| 50 | +#include "Framework/ASoA.h" |
| 51 | +#include "Framework/AnalysisDataModel.h" |
| 52 | +#include "Framework/AnalysisTask.h" |
| 53 | +#include "Framework/O2DatabasePDGPlugin.h" |
| 54 | +#include "Framework/HistogramRegistry.h" |
| 55 | + |
| 56 | +#include "Common/Core/TrackSelection.h" |
| 57 | +#include "Common/Core/TrackSelectionDefaults.h" |
| 58 | + |
| 59 | +#include "Common/DataModel/EventSelection.h" |
| 60 | +#include "Common/DataModel/TrackSelectionTables.h" |
| 61 | + |
| 62 | +#include "PWGJE/Core/FastJetUtilities.h" |
| 63 | +#include "PWGJE/Core/JetFinder.h" |
| 64 | +#include "PWGJE/Core/JetFindingUtilities.h" |
| 65 | +#include "PWGJE/DataModel/Jet.h" |
| 66 | + |
| 67 | +#include "PWGJE/Core/JetDerivedDataUtilities.h" |
| 68 | +#include "EventFiltering/filterTables.h" |
| 69 | + |
| 70 | +using namespace o2; |
| 71 | +using namespace o2::framework; |
| 72 | +using namespace o2::framework::expressions; |
| 73 | + |
| 74 | +namespace o2::aod |
| 75 | +{ |
| 76 | +DECLARE_SOA_TABLE(MyCollisions, "AOD", "MYCOLLISION", //! vertex information of collision |
| 77 | + o2::soa::Index<>, collision::PosZ); |
| 78 | +using MyCollision = MyCollisions::iterator; |
| 79 | + |
| 80 | + namespace myTable |
| 81 | + { |
| 82 | + DECLARE_SOA_INDEX_COLUMN(MyCollision, mycollision); |
| 83 | + |
| 84 | + DECLARE_SOA_COLUMN(Psi2, psi2, Float_t); |
| 85 | + DECLARE_SOA_COLUMN(Psi3, psi3, Float_t); |
| 86 | + DECLARE_SOA_COLUMN(EvtPlRes, evtplres, Float_t); |
| 87 | + } |
| 88 | + DECLARE_SOA_TABLE(MyTable, "AOD", "MYTABLE", o2::soa::Index<>, |
| 89 | + myTable::Psi2, |
| 90 | + myTable::Psi3, |
| 91 | + myTable::EvtPlRes); |
| 92 | +} |
| 93 | + |
| 94 | +///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// |
| 95 | +struct jetChargedV2EPTable { |
| 96 | + Produces<o2::aod::MyTable> myTable; |
| 97 | + Produces<o2::aod::MyCollision> outputCollisions; |
| 98 | + |
| 99 | + HistogramRegistry registry; |
| 100 | + |
| 101 | + |
| 102 | + Configurable<std::string> eventSelections{"eventSelections", "sel8", "choose event selection"}; |
| 103 | + Configurable<std::string> trackSelections{"trackSelections", "globalTracks", "set track selections"}; |
| 104 | + |
| 105 | + Configurable<float> vertexZCut{"vertexZCut", 10.0f, "Accepted z-vertex range"}; |
| 106 | + Configurable<float> centralityMin{"centralityMin", -999.0, "minimum centrality"}; |
| 107 | + Configurable<float> centralityMax{"centralityMax", 999.0, "maximum centrality"}; |
| 108 | + Configurable<float> trackPtMin{"trackPtMin", 0.15, "minimum pT acceptance for tracks"}; |
| 109 | + Configurable<float> trackPtMax{"trackPtMax", 1000., "maximum pT acceptance for tracks"}; |
| 110 | + Configurable<float> trackEtaMin{"trackEtaMin", -0.9, "minimum eta acceptance for tracks"}; |
| 111 | + Configurable<float> trackEtaMax{"trackEtaMax", 0.9, "maximum eta acceptance for tracks"}; |
| 112 | + |
| 113 | + Configurable<float> jetAreaFractionMin{"jetAreaFractionMin", -99.0, "used to make a cut on the jet areas"}; |
| 114 | + Configurable<float> leadingConstituentPtMin{"leadingConstituentPtMin", -99.0, "minimum pT selection on jet constituent"}; |
| 115 | + Configurable<float> jetPtMin{"jetPtMin", 0.15, "minimum pT acceptance for jets"}; |
| 116 | + Configurable<float> jetPtMax{"jetPtMax", 200.0, "maximum pT acceptance for jets"}; |
| 117 | + |
| 118 | + Configurable<int> trackOccupancyInTimeRangeMax{"trackOccupancyInTimeRangeMax", 999999, "maximum occupancy of tracks in neighbouring collisions in a given time range; only applied to reconstructed collisions (data and mcd jets), not mc collisions (mcp jets)"}; |
| 119 | + Configurable<int> trackOccupancyInTimeRangeMin{"trackOccupancyInTimeRangeMin", -999999, "minimum occupancy of tracks in neighbouring collisions in a given time range; only applied to reconstructed collisions (data and mcd jets), not mc collisions (mcp jets)"}; |
| 120 | + |
| 121 | + //=====================< evt pln >=====================// |
| 122 | + Configurable<bool> cfgAddEvtSel{"cfgAddEvtSel", true, "event selection"}; |
| 123 | + Configurable<std::vector<int>> cfgnMods{"cfgnMods", {2}, "Modulation of interest"}; |
| 124 | + Configurable<int> cfgnTotalSystem{"cfgnTotalSystem", 7, "total qvect 56 or number"}; |
| 125 | + Configurable<std::string> cfgDetName{"cfgDetName", "FT0M", "The name of detector to be analyzed"}; |
| 126 | + Configurable<std::string> cfgRefAName{"cfgRefAName", "TPCpos", "The name of detector for reference A"}; |
| 127 | + Configurable<std::string> cfgRefBName{"cfgRefBName", "TPCneg", "The name of detector for reference B"}; |
| 128 | + |
| 129 | + EventPlaneHelper helperEP; |
| 130 | + int DetId; |
| 131 | + int RefAId; |
| 132 | + int RefBId; |
| 133 | + |
| 134 | + template <typename T> |
| 135 | + int GetDetId(const T& name) |
| 136 | + { |
| 137 | + if (name.value == "BPos" || name.value == "BNeg" || name.value == "BTot") { |
| 138 | + LOGF(warning, "Using deprecated label: %s. Please use TPCpos, TPCneg, TPCall instead.", name.value); |
| 139 | + } |
| 140 | + if (name.value == "FT0C") { |
| 141 | + return 0; |
| 142 | + } else if (name.value == "FT0A") { |
| 143 | + return 1; |
| 144 | + } else if (name.value == "FT0M") { |
| 145 | + return 2; |
| 146 | + } else if (name.value == "FV0A") { |
| 147 | + return 3; |
| 148 | + } else if (name.value == "TPCpos" || name.value == "BPos") { |
| 149 | + return 4; |
| 150 | + } else if (name.value == "TPCneg" || name.value == "BNeg") { |
| 151 | + return 5; |
| 152 | + } else if (name.value == "TPCall" || name.value == "BTot") { |
| 153 | + return 6; |
| 154 | + } else { |
| 155 | + return 0; |
| 156 | + } |
| 157 | + } |
| 158 | + //=====================< evt pln | end >=====================// |
| 159 | + |
| 160 | + Configurable<float> selectedJetsRadius{"selectedJetsRadius", 0.4, "resolution parameter for histograms without radius"}; |
| 161 | + |
| 162 | + std::vector<double> jetPtBins; |
| 163 | + std::vector<double> jetPtBinsRhoAreaSub; |
| 164 | + |
| 165 | + int eventSelectiont = -1; |
| 166 | + int trackSelectiont = -1; |
| 167 | + |
| 168 | + void init(o2::framework::InitContext&) |
| 169 | + { |
| 170 | + DetId = GetDetId(cfgDetName); |
| 171 | + RefAId = GetDetId(cfgRefAName); |
| 172 | + RefBId = GetDetId(cfgRefBName); |
| 173 | + if (DetId == RefAId || DetId == RefBId || RefAId == RefBId) { |
| 174 | + LOGF(info, "Wrong detector configuration \n The FT0C will be used to get Q-Vector \n The TPCpos and TPCneg will be used as reference systems"); |
| 175 | + DetId = 0; |
| 176 | + RefAId = 4; |
| 177 | + RefBId = 5; |
| 178 | + } |
| 179 | + eventSelectiont = jetderiveddatautilities::initialiseEventSelection(static_cast<std::string>(eventSelections)); |
| 180 | + trackSelectiont = jetderiveddatautilities::initialiseTrackSelection(static_cast<std::string>(trackSelections)); |
| 181 | + } |
| 182 | + |
| 183 | + Filter trackCuts = (aod::jtrack::pt >= trackPtMin && aod::jtrack::pt < trackPtMax && aod::jtrack::eta > trackEtaMin && aod::jtrack::eta < trackEtaMax); |
| 184 | + Filter eventCuts = (nabs(aod::jcollision::posZ) < vertexZCut && aod::jcollision::centrality >= centralityMin && aod::jcollision::centrality < centralityMax); |
| 185 | + |
| 186 | + void processJetV2Table(soa::Filtered<soa::Join<aod::JetCollisions, aod::BkgChargedRhos, aod::Qvectors>>::iterator const& collision, |
| 187 | + soa::Join<aod::ChargedJets, aod::ChargedJetConstituents> const& jets, |
| 188 | + aod::JetTracks const& tracks) |
| 189 | + { |
| 190 | + if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) { |
| 191 | + return; |
| 192 | + } |
| 193 | + if (!jetderiveddatautilities::selectCollision(collision, eventSelectiont)) { |
| 194 | + return; |
| 195 | + } |
| 196 | + outputCollisions(collision.posZ()); |
| 197 | + //=====================< evt pln [n=2->\Psi_2, n=3->\Psi_3] >=====================// |
| 198 | + double ep2 = 0; |
| 199 | + double ep3 = 0; |
| 200 | + double res2 = 0; |
| 201 | + |
| 202 | + for (uint i = 0; i < cfgnMods->size(); i++) { |
| 203 | + int nmode = cfgnMods->at(i); |
| 204 | + int DetInd = DetId * 4 + cfgnTotalSystem * 4 * (nmode - 2); |
| 205 | + if (nmode == 2) { |
| 206 | + if (collision.qvecAmp()[DetId] > 1e-8) { |
| 207 | + ep2 = helperEP.GetEventPlane(collision.qvecRe()[DetInd + 3], collision.qvecIm()[DetInd + 3], nmode); |
| 208 | + res2 = helperEP.GetResolution(helperEP.GetEventPlane(collision.qvecRe()[DetInd + 3], collision.qvecIm()[DetInd + 3], nmode), helperEP.GetEventPlane(collision.qvecRe()[DetInd + 3], collision.qvecIm()[DetInd + 3], nmode), nmode); |
| 209 | + } |
| 210 | + } else if (nmode == 3) { |
| 211 | + if (collision.qvecAmp()[DetId] > 1e-8) { |
| 212 | + ep3 = helperEP.GetEventPlane(collision.qvecRe()[DetInd + 3], collision.qvecIm()[DetInd + 3], nmode); |
| 213 | + } |
| 214 | + } |
| 215 | + } |
| 216 | + myTable(ep2, ep3, res2); |
| 217 | + LOGF(info, "table producer Psi2 as %f, Psi3 as %f, Res2 as %f", ep2, ep3, res2); |
| 218 | + } |
| 219 | + PROCESS_SWITCH(jetChargedV2EPTable, processJetV2Table, "Jet V2 table for fit", true); |
| 220 | +}; |
| 221 | + |
| 222 | +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) |
| 223 | +{ |
| 224 | + return WorkflowSpec{adaptAnalysisTask<jetChargedV2EPTable>(cfgc, TaskName{"jet-charged-v2ept"})}; |
| 225 | +} |
0 commit comments