|
| 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 | +#include <vector> |
| 13 | +#include <TStopwatch.h> |
| 14 | +#include "DataFormatsGlobalTracking/RecoContainer.h" |
| 15 | +#include "DataFormatsGlobalTracking/RecoContainerCreateTracksVariadic.h" |
| 16 | +#include "DetectorsVertexing/SVertexerParams.h" |
| 17 | +#include "ReconstructionDataFormats/TrackTPCITS.h" |
| 18 | +#include "ReconstructionDataFormats/V0.h" |
| 19 | +#include "ReconstructionDataFormats/GlobalTrackID.h" |
| 20 | +#include "DataFormatsCalibration/MeanVertexObject.h" |
| 21 | +#include "DetectorsBase/Propagator.h" |
| 22 | +#include "DetectorsBase/GeometryManager.h" |
| 23 | +#include "SimulationDataFormat/MCEventLabel.h" |
| 24 | +#include "SimulationDataFormat/MCUtils.h" |
| 25 | +#include "CommonDataFormat/BunchFilling.h" |
| 26 | +#include "CommonUtils/NameConf.h" |
| 27 | +#include "DataFormatsFT0/RecPoints.h" |
| 28 | +#include "Framework/ConfigParamRegistry.h" |
| 29 | +#include "Framework/CCDBParamSpec.h" |
| 30 | +#include "FT0Reconstruction/InteractionTag.h" |
| 31 | +#include "ITSMFTBase/DPLAlpideParam.h" |
| 32 | +#include "DetectorsCommonDataFormats/DetID.h" |
| 33 | +#include "DetectorsBase/GRPGeomHelper.h" |
| 34 | +#include "GlobalTrackingStudy/TrackingStudy.h" |
| 35 | +#include "TPCBase/ParameterElectronics.h" |
| 36 | +#include "ReconstructionDataFormats/PrimaryVertex.h" |
| 37 | +#include "DataFormatsFT0/RecPoints.h" |
| 38 | +#include "DataFormatsITSMFT/TrkClusRef.h" |
| 39 | +#include "CommonUtils/TreeStreamRedirector.h" |
| 40 | +#include "ReconstructionDataFormats/VtxTrackRef.h" |
| 41 | +#include "ReconstructionDataFormats/DCA.h" |
| 42 | +#include "Steer/MCKinematicsReader.h" |
| 43 | +#include "DCAFitter/DCAFitterN.h" |
| 44 | +#include "MathUtils/fit.h" |
| 45 | + |
| 46 | +namespace o2::svstudy |
| 47 | +{ |
| 48 | + |
| 49 | +using namespace o2::framework; |
| 50 | +using DetID = o2::detectors::DetID; |
| 51 | +using DataRequest = o2::globaltracking::DataRequest; |
| 52 | + |
| 53 | +using PVertex = o2::dataformats::PrimaryVertex; |
| 54 | +using V2TRef = o2::dataformats::VtxTrackRef; |
| 55 | +using VTIndex = o2::dataformats::VtxTrackIndex; |
| 56 | +using GTrackID = o2::dataformats::GlobalTrackID; |
| 57 | +using TBracket = o2::math_utils::Bracketf_t; |
| 58 | +using V0ID = o2::dataformats::V0Index; |
| 59 | + |
| 60 | +using timeEst = o2::dataformats::TimeStampWithError<float, float>; |
| 61 | + |
| 62 | +class SVStudySpec : public Task |
| 63 | +{ |
| 64 | + public: |
| 65 | + SVStudySpec(std::shared_ptr<DataRequest> dr, std::shared_ptr<o2::base::GRPGeomRequest> gr, GTrackID::mask_t src, bool useMC) |
| 66 | + : mDataRequest(dr), mGGCCDBRequest(gr), mTracksSrc(src), mUseMC(useMC) {} |
| 67 | + ~SVStudySpec() final = default; |
| 68 | + void init(InitContext& ic) final; |
| 69 | + void run(ProcessingContext& pc) final; |
| 70 | + void endOfStream(EndOfStreamContext& ec) final; |
| 71 | + void finaliseCCDB(ConcreteDataMatcher& matcher, void* obj) final; |
| 72 | + void process(o2::globaltracking::RecoContainer& recoData); |
| 73 | + |
| 74 | + private: |
| 75 | + void updateTimeDependentParams(ProcessingContext& pc); |
| 76 | + bool refitV0(const V0ID& id, o2::dataformats::V0& v0, o2::globaltracking::RecoContainer& recoData); |
| 77 | + std::shared_ptr<DataRequest> mDataRequest; |
| 78 | + std::shared_ptr<o2::base::GRPGeomRequest> mGGCCDBRequest; |
| 79 | + bool mUseMC{false}; ///< MC flag |
| 80 | + std::unique_ptr<o2::utils::TreeStreamRedirector> mDBGOut; |
| 81 | + bool mSelK0 = false; |
| 82 | + bool mRefit = false; |
| 83 | + float mMaxEta = 0.8; |
| 84 | + float mBz = 0; |
| 85 | + GTrackID::mask_t mTracksSrc{}; |
| 86 | + o2::vertexing::DCAFitterN<2> mFitterV0; |
| 87 | + o2::steer::MCKinematicsReader mcReader; // reader of MC information |
| 88 | +}; |
| 89 | + |
| 90 | +void SVStudySpec::init(InitContext& ic) |
| 91 | +{ |
| 92 | + o2::base::GRPGeomHelper::instance().setRequest(mGGCCDBRequest); |
| 93 | + mDBGOut = std::make_unique<o2::utils::TreeStreamRedirector>("svStudy.root", "recreate"); |
| 94 | + mRefit = ic.options().get<bool>("refit"); |
| 95 | + mSelK0 = ic.options().get<bool>("sel-k0"); |
| 96 | + mMaxEta = ic.options().get<float>("max-eta"); |
| 97 | +} |
| 98 | + |
| 99 | +void SVStudySpec::run(ProcessingContext& pc) |
| 100 | +{ |
| 101 | + o2::globaltracking::RecoContainer recoData; |
| 102 | + recoData.collectData(pc, *mDataRequest.get()); // select tracks of needed type, with minimal cuts, the real selected will be done in the vertexer |
| 103 | + updateTimeDependentParams(pc); // Make sure this is called after recoData.collectData, which may load some conditions |
| 104 | + process(recoData); |
| 105 | +} |
| 106 | + |
| 107 | +void SVStudySpec::updateTimeDependentParams(ProcessingContext& pc) |
| 108 | +{ |
| 109 | + o2::base::GRPGeomHelper::instance().checkUpdates(pc); |
| 110 | + static bool initOnceDone = false; |
| 111 | + if (!initOnceDone) { // this params need to be queried only once |
| 112 | + initOnceDone = true; |
| 113 | + const auto& svparam = o2::vertexing::SVertexerParams::Instance(); |
| 114 | + // Note: reading of the ITS AlpideParam needed for ITS timing is done by the RecoContainer |
| 115 | + mFitterV0.setBz(mBz); |
| 116 | + mFitterV0.setUseAbsDCA(svparam.useAbsDCA); |
| 117 | + mFitterV0.setPropagateToPCA(false); |
| 118 | + mFitterV0.setMaxR(svparam.maxRIni); |
| 119 | + mFitterV0.setMinParamChange(svparam.minParamChange); |
| 120 | + mFitterV0.setMinRelChi2Change(svparam.minRelChi2Change); |
| 121 | + mFitterV0.setMaxDZIni(svparam.maxDZIni); |
| 122 | + mFitterV0.setMaxDXYIni(svparam.maxDXYIni); |
| 123 | + mFitterV0.setMaxChi2(svparam.maxChi2); |
| 124 | + mFitterV0.setMatCorrType(o2::base::Propagator::MatCorrType(svparam.matCorr)); |
| 125 | + mFitterV0.setUsePropagator(svparam.usePropagator); |
| 126 | + mFitterV0.setRefitWithMatCorr(svparam.refitWithMatCorr); |
| 127 | + mFitterV0.setMaxStep(svparam.maxStep); |
| 128 | + mFitterV0.setMaxSnp(svparam.maxSnp); |
| 129 | + mFitterV0.setMinXSeed(svparam.minXSeed); |
| 130 | + } |
| 131 | + mBz = o2::base::Propagator::Instance()->getNominalBz(); |
| 132 | + mFitterV0.setBz(mBz); |
| 133 | +} |
| 134 | + |
| 135 | +void SVStudySpec::process(o2::globaltracking::RecoContainer& recoData) |
| 136 | +{ |
| 137 | + auto v0s = recoData.getV0s(); |
| 138 | + auto v0IDs = recoData.getV0sIdx(); |
| 139 | + bool refit = mRefit || (v0s.size() < v0IDs.size()); |
| 140 | + int nv0 = v0IDs.size(); |
| 141 | + o2::dataformats::V0 v0ref; |
| 142 | + o2::track::TrackParCov dummyTr{}; |
| 143 | + const o2::track::TrackParCov* tpcTracks[2] = {&dummyTr, &dummyTr}; |
| 144 | + int nclTPC[2] = {0, 0}, nclITS[2] = {0, 0}, itsPatt[2] = {0, 0}; |
| 145 | + static int tfID = 0; |
| 146 | + |
| 147 | + for (int iv = 0; iv < nv0; iv++) { |
| 148 | + const auto& v0id = v0IDs[iv]; |
| 149 | + if (mRefit && !refitV0(v0id, v0ref, recoData)) { |
| 150 | + continue; |
| 151 | + } |
| 152 | + const auto& v0 = mRefit ? v0ref : v0s[iv]; |
| 153 | + if (mMaxEta > std::abs(v0.getEta())) { |
| 154 | + continue; |
| 155 | + } |
| 156 | + if (mSelK0 && std::abs(std::sqrt(v0.calcMass2AsK0()) - 0.497) > 0.1) { |
| 157 | + continue; |
| 158 | + } |
| 159 | + for (int ip = 0; ip < 2; ip++) { |
| 160 | + auto gid = v0id.getProngID(ip); |
| 161 | + |
| 162 | + // get TPC tracks, if any |
| 163 | + tpcTracks[ip] = &dummyTr; |
| 164 | + nclTPC[ip] = 0; |
| 165 | + if (gid.includesDet(DetID::TPC)) { |
| 166 | + const auto& tpcTr = recoData.getTPCTrack(recoData.getTPCContributorGID(gid)); |
| 167 | + tpcTracks[ip] = &tpcTr; |
| 168 | + nclTPC[ip] = tpcTr.getNClusters(); |
| 169 | + } |
| 170 | + // get ITS tracks, if any |
| 171 | + nclITS[ip] = itsPatt[ip] = 0; |
| 172 | + if (gid.includesDet(DetID::ITS)) { |
| 173 | + auto gidITS = recoData.getITSContributorGID(gid); |
| 174 | + if (gidITS.getSource() == GTrackID::ITS) { |
| 175 | + const auto& itsTr = recoData.getITSTrack(recoData.getITSContributorGID(gid)); |
| 176 | + nclITS[ip] = itsTr.getNClusters(); |
| 177 | + for (int il = 0; il < 7; il++) { |
| 178 | + if (itsTr.hasHitOnLayer(il)) { |
| 179 | + itsPatt[ip] |= 0x1 << il; |
| 180 | + } |
| 181 | + } |
| 182 | + } else { |
| 183 | + const auto& itsTrf = recoData.getITSABRefs()[gidITS]; |
| 184 | + nclITS[ip] = itsTrf.getNClusters(); |
| 185 | + for (int il = 0; il < 7; il++) { |
| 186 | + if (itsTrf.hasHitOnLayer(il)) { |
| 187 | + itsPatt[ip] |= 0x1 << il; |
| 188 | + } |
| 189 | + } |
| 190 | + } |
| 191 | + } |
| 192 | + } |
| 193 | + (*mDBGOut) << "tfinfo" |
| 194 | + << "orbit=" << recoData.startIR.orbit << "tfID=" << tfID << "\n"; |
| 195 | + (*mDBGOut) << "v0s" |
| 196 | + << "v0=" << v0 << "v0ID=" << v0id << "tpc0=" << *tpcTracks[0] << "tpc1=" << *tpcTracks[1] |
| 197 | + << "nclTPC0=" << nclTPC[0] << "nclTPC1=" << nclTPC[1] << "nclITS0=" << nclITS[0] << "nclITS1=" << nclITS[1] << "itsPatt0=" << itsPatt[0] << "itsPatt1=" << itsPatt[1] << "\n"; |
| 198 | + } |
| 199 | + tfID++; |
| 200 | +} |
| 201 | + |
| 202 | +bool SVStudySpec::refitV0(const V0ID& id, o2::dataformats::V0& v0, o2::globaltracking::RecoContainer& recoData) |
| 203 | +{ |
| 204 | + auto seedP = recoData.getTrackParam(id.getProngID(0)); |
| 205 | + auto seedN = recoData.getTrackParam(id.getProngID(1)); |
| 206 | + bool isTPConly = (id.getProngID(0).getSource() == GTrackID::TPC) || (id.getProngID(1).getSource() == GTrackID::TPC); |
| 207 | + const auto& svparam = o2::vertexing::SVertexerParams::Instance(); |
| 208 | + if (svparam.mTPCTrackPhotonTune && isTPConly) { |
| 209 | + mFitterV0.setMaxDZIni(svparam.mTPCTrackMaxDZIni); |
| 210 | + mFitterV0.setMaxDXYIni(svparam.mTPCTrackMaxDXYIni); |
| 211 | + mFitterV0.setMaxChi2(svparam.mTPCTrackMaxChi2); |
| 212 | + mFitterV0.setCollinear(true); |
| 213 | + } |
| 214 | + int nCand = mFitterV0.process(seedP, seedN); |
| 215 | + if (svparam.mTPCTrackPhotonTune && isTPConly) { // restore |
| 216 | + // Reset immediately to the defaults |
| 217 | + mFitterV0.setMaxDZIni(svparam.maxDZIni); |
| 218 | + mFitterV0.setMaxDXYIni(svparam.maxDXYIni); |
| 219 | + mFitterV0.setMaxChi2(svparam.maxChi2); |
| 220 | + mFitterV0.setCollinear(false); |
| 221 | + } |
| 222 | + if (nCand == 0) { // discard this pair |
| 223 | + return false; |
| 224 | + } |
| 225 | + const int cand = 0; |
| 226 | + if (!mFitterV0.isPropagateTracksToVertexDone(cand) && !mFitterV0.propagateTracksToVertex(cand)) { |
| 227 | + return false; |
| 228 | + } |
| 229 | + const auto& trPProp = mFitterV0.getTrack(0, cand); |
| 230 | + const auto& trNProp = mFitterV0.getTrack(1, cand); |
| 231 | + std::array<float, 3> pP{}, pN{}; |
| 232 | + trPProp.getPxPyPzGlo(pP); |
| 233 | + trNProp.getPxPyPzGlo(pN); |
| 234 | + std::array<float, 3> pV0 = {pP[0] + pN[0], pP[1] + pN[1], pP[2] + pN[2]}; |
| 235 | + auto p2V0 = pV0[0] * pV0[0] + pV0[1] * pV0[1] + pV0[2] * pV0[2]; |
| 236 | + const auto& pv = recoData.getPrimaryVertex(id.getVertexID()); |
| 237 | + const auto v0XYZ = mFitterV0.getPCACandidatePos(cand); |
| 238 | + float dx = v0XYZ[0] - pv.getX(), dy = v0XYZ[1] - pv.getY(), dz = v0XYZ[2] - pv.getZ(), prodXYZv0 = dx * pV0[0] + dy * pV0[1] + dz * pV0[2]; |
| 239 | + float cosPA = prodXYZv0 / std::sqrt((dx * dx + dy * dy + dz * dz) * p2V0); |
| 240 | + new (&v0) o2::dataformats::V0(v0XYZ, pV0, mFitterV0.calcPCACovMatrixFlat(cand), trPProp, trNProp); |
| 241 | + v0.setDCA(mFitterV0.getChi2AtPCACandidate(cand)); |
| 242 | + v0.setCosPA(cosPA); |
| 243 | + return true; |
| 244 | +} |
| 245 | + |
| 246 | +void SVStudySpec::endOfStream(EndOfStreamContext& ec) |
| 247 | +{ |
| 248 | + mDBGOut.reset(); |
| 249 | +} |
| 250 | + |
| 251 | +void SVStudySpec::finaliseCCDB(ConcreteDataMatcher& matcher, void* obj) |
| 252 | +{ |
| 253 | + if (o2::base::GRPGeomHelper::instance().finaliseCCDB(matcher, obj)) { |
| 254 | + return; |
| 255 | + } |
| 256 | +} |
| 257 | + |
| 258 | +DataProcessorSpec getSVStudySpec(GTrackID::mask_t srcTracks, bool useMC) |
| 259 | +{ |
| 260 | + std::vector<OutputSpec> outputs; |
| 261 | + auto dataRequest = std::make_shared<DataRequest>(); |
| 262 | + |
| 263 | + dataRequest->requestTracks(srcTracks, useMC); |
| 264 | + dataRequest->requestPrimaryVertertices(useMC); |
| 265 | + dataRequest->requestSecondaryVertices(useMC); |
| 266 | + dataRequest->inputs.emplace_back("meanvtx", "GLO", "MEANVERTEX", 0, Lifetime::Condition, ccdbParamSpec("GLO/Calib/MeanVertex", {}, 1)); |
| 267 | + auto ggRequest = std::make_shared<o2::base::GRPGeomRequest>(false, // orbitResetTime |
| 268 | + false, // GRPECS=true |
| 269 | + false, // GRPLHCIF |
| 270 | + true, // GRPMagField |
| 271 | + true, // askMatLUT |
| 272 | + o2::base::GRPGeomRequest::None, // geometry |
| 273 | + dataRequest->inputs, |
| 274 | + true); |
| 275 | + |
| 276 | + return DataProcessorSpec{ |
| 277 | + "sv-study", |
| 278 | + dataRequest->inputs, |
| 279 | + outputs, |
| 280 | + AlgorithmSpec{adaptFromTask<SVStudySpec>(dataRequest, ggRequest, srcTracks, useMC)}, |
| 281 | + Options{ |
| 282 | + {"refit", VariantType::Bool, false, {"refit SVertices"}}, |
| 283 | + {"sel-k0", VariantType::Bool, false, {"select K0s only"}}, |
| 284 | + {"max-eta", VariantType::Float, 1.2f, {"Cut on track eta"}}, |
| 285 | + }}; |
| 286 | +} |
| 287 | + |
| 288 | +} // namespace o2::svstudy |
0 commit comments