Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:13:03

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008 #pragma once
0009 
0010 #include "Acts/Definitions/Algebra.hpp"
0011 #include "Acts/Definitions/Units.hpp"
0012 #include "Acts/Surfaces/PerigeeSurface.hpp"
0013 #include "Acts/Tests/CommonHelpers/DataDirectory.hpp"
0014 #include "Acts/Vertexing/Vertex.hpp"
0015 
0016 #include <fstream>
0017 #include <iterator>
0018 #include <regex>
0019 
0020 namespace Acts::Test {
0021 
0022 using namespace Acts::UnitLiterals;
0023 using Covariance = BoundSquareMatrix;
0024 
0025 enum VertexCsvData { BeamSpotData, VerticesData, TracksData };
0026 
0027 /// @brief Helper struct to store reference vertex related information
0028 struct VertexInfo {
0029   // The position
0030   Vector3 position;
0031   // The covariance
0032   SquareMatrix3 covariance;
0033   // Number of tracks
0034   int nTracks = 0;
0035   // Weight of first track
0036   double trk1Weight = 0;
0037   // Vertex compatibility value of first track
0038   double trk1Comp = 0;
0039   // Chi2 of first track
0040   double trk1Chi2 = 0;
0041 };
0042 
0043 inline std::tuple<Vertex, std::vector<VertexInfo>,
0044                   std::vector<BoundTrackParameters>>
0045 readTracksAndVertexCSV(const std::string& toolString,
0046                        const std::string& fileBase = "vertexing_event_mu20") {
0047   const auto beamspotDataPath =
0048       Acts::Test::getDataPath(fileBase + "_beamspot.csv");
0049   const auto tracksDataPath = Acts::Test::getDataPath(fileBase + "_tracks.csv");
0050   const auto verticesDataPath =
0051       Acts::Test::getDataPath(fileBase + "_vertices_" + toolString + ".csv");
0052 
0053   const std::regex comma(",");
0054 
0055   // Open source files
0056   std::ifstream beamspotData(beamspotDataPath);
0057   std::ifstream tracksData(tracksDataPath);
0058   std::ifstream verticesData(verticesDataPath);
0059 
0060   // String to store the read lines
0061   std::string line{};
0062 
0063   std::shared_ptr<PerigeeSurface> perigeeSurface;
0064   std::vector<BoundTrackParameters> tracks;
0065   std::vector<VertexInfo> vertices;
0066   Vertex beamspotConstraint;
0067 
0068   // Read in beamspot data
0069   std::getline(beamspotData, line);  // skip header
0070   while (beamspotData && std::getline(beamspotData, line)) {
0071     // Tokenize line and store result in vector
0072     std::vector<std::string> row{
0073         std::sregex_token_iterator(line.begin(), line.end(), comma, -1),
0074         std::sregex_token_iterator()};
0075 
0076     Vector3 beamspotPos;
0077     SquareMatrix3 beamspotCov;
0078     beamspotPos << std::stod(row[0]) * (1_mm), std::stod(row[1]) * (1_mm),
0079         std::stod(row[2]) * (1_mm);
0080     beamspotCov << std::stod(row[3]), 0, 0, 0, std::stod(row[4]), 0, 0, 0,
0081         std::stod(row[5]);
0082     beamspotConstraint.setPosition(beamspotPos);
0083     beamspotConstraint.setCovariance(beamspotCov);
0084     perigeeSurface = Surface::makeShared<PerigeeSurface>(beamspotPos);
0085   }
0086 
0087   // Read in track data
0088   std::getline(tracksData, line);  // skip header
0089   while (tracksData && std::getline(tracksData, line)) {
0090     // Tokenize line and store result in vector
0091     std::vector<std::string> row{
0092         std::sregex_token_iterator(line.begin(), line.end(), comma, -1),
0093         std::sregex_token_iterator()};
0094 
0095     BoundVector params;
0096     params << std::stod(row[0]), std::stod(row[1]), std::stod(row[2]),
0097         std::stod(row[3]), std::stod(row[4]) * 1. / (1_MeV), std::stod(row[5]);
0098     Covariance covMat;
0099     covMat << std::stod(row[6]), std::stod(row[7]), std::stod(row[8]),
0100         std::stod(row[9]), std::stod(row[10]) * 1. / (1_MeV),
0101         std::stod(row[11]), std::stod(row[7]), std::stod(row[12]),
0102         std::stod(row[13]), std::stod(row[14]),
0103         std::stod(row[15]) * 1. / (1_MeV), std::stod(row[16]),
0104         std::stod(row[8]), std::stod(row[13]), std::stod(row[17]),
0105         std::stod(row[18]), std::stod(row[19]) * 1. / (1_MeV),
0106         std::stod(row[20]), std::stod(row[9]), std::stod(row[14]),
0107         std::stod(row[18]), std::stod(row[21]),
0108         std::stod(row[22]) * 1. / (1_MeV), std::stod(row[23]),
0109         std::stod(row[10]) * 1. / (1_MeV), std::stod(row[15]) * 1. / (1_MeV),
0110         std::stod(row[19]) * 1. / (1_MeV), std::stod(row[22]) * 1. / (1_MeV),
0111         std::stod(row[24]) * 1. / (1_MeV * 1_MeV),
0112         std::stod(row[25]) * 1. / (1_MeV), std::stod(row[11]),
0113         std::stod(row[16]), std::stod(row[20]), std::stod(row[23]),
0114         std::stod(row[25]) * 1. / (1_MeV), std::stod(row[26]);
0115 
0116     // TODO we do not have a hypothesis at hand here. defaulting to pion
0117     tracks.emplace_back(perigeeSurface, params, std::move(covMat),
0118                         ParticleHypothesis::pion());
0119   }
0120 
0121   // Read in reference vertex data
0122   std::getline(verticesData, line);  // skip header
0123   while (verticesData && std::getline(verticesData, line)) {
0124     // Tokenize line and store result in vector
0125     std::vector<std::string> row{
0126         std::sregex_token_iterator(line.begin(), line.end(), comma, -1),
0127         std::sregex_token_iterator()};
0128 
0129     Vector3 pos;
0130     pos << std::stod(row[0]) * (1_mm), std::stod(row[1]) * (1_mm),
0131         std::stod(row[2]) * (1_mm);
0132     SquareMatrix3 cov;
0133     cov << std::stod(row[3]), std::stod(row[4]), std::stod(row[5]),
0134         std::stod(row[6]), std::stod(row[7]), std::stod(row[8]),
0135         std::stod(row[9]), std::stod(row[10]), std::stod(row[11]);
0136     VertexInfo vertexInfo;
0137     vertexInfo.position = pos;
0138     vertexInfo.covariance = cov;
0139     vertexInfo.nTracks = std::stoi(row[12]);
0140     vertexInfo.trk1Weight = std::stod(row[13]);
0141     vertexInfo.trk1Comp = std::stod(row[14]);
0142     vertexInfo.trk1Chi2 = std::stod(row[15]);
0143     vertices.push_back(vertexInfo);
0144   }
0145 
0146   return {beamspotConstraint, vertices, tracks};
0147 }
0148 
0149 }  // namespace Acts::Test