src/MagneticField/TextMagneticFieldIo.cpp

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
0009 #include "Acts/MagneticField/TextMagneticFieldIo.hpp"
0010
0011 #include "Acts/MagneticField/BFieldMapUtils.hpp"
0012
0013 #include <fstream>
0014 #include <iostream>
0015 #include <vector>
0016
0017 #include <boost/algorithm/string.hpp>
0018
0019 namespace {
0020 constexpr std::size_t kDefaultSize = 1 << 15;
0021
0022 bool ignoreLine(const std::string& line) {
0023   return line.empty() || line[0] == '%' || line[0] == '#' ||
0024          std::isalpha(line[0]) != 0 ||
0025          line.find_first_not_of(' ') == std::string::npos;
0026 }
0027
0028 }  // namespace
0029
0030 Acts::InterpolatedBFieldMap<
0031     Acts::Grid<Acts::Vector2, Acts::Axis<Acts::AxisType::Equidistant>,
0032                Acts::Axis<Acts::AxisType::Equidistant>>>
0033 Acts::makeMagneticFieldMapRzFromText(
0034     const std::function<std::size_t(std::array<std::size_t, 2> binsRZ,
0035                                     std::array<std::size_t, 2> nBinsRZ)>&
0036         localToGlobalBin,
0037     const std::string& fieldMapFile, double lengthUnit, double BFieldUnit,
0038     bool firstQuadrant, const std::string& delimiter) {
0039   /// [1] Read in field map file
0040   // Grid position points in r and z
0041   std::vector<double> rPos;
0042   std::vector<double> zPos;
0043   // components of magnetic field on grid points
0044   std::vector<Vector2> bField;
0045   // reserve estimated size
0046   rPos.reserve(kDefaultSize);
0047   zPos.reserve(kDefaultSize);
0048   bField.reserve(kDefaultSize);
0049   // [1] Read in file and fill values
0050   std::ifstream map_file(fieldMapFile.c_str(), std::ios::in);
0051   std::string line;
0052   double r = 0., z = 0.;
0053   double br = 0., bz = 0.;
0054   while (std::getline(map_file, line)) {
0055     if (ignoreLine(line)) {
0056       continue;
0057     }
0058
0059     std::istringstream tmp(line);
0060     if (delimiter.empty()) {
0061       tmp >> r >> z >> br >> bz;
0062     } else {
0063       std::vector<std::string> tokens;
0064       boost::split(tokens, line, boost::is_any_of(delimiter));
0065       if (tokens.size() == 4) {
0066         r = std::stod(tokens[0]);
0067         z = std::stod(tokens[1]);
0068         br = std::stod(tokens[2]);
0069         bz = std::stod(tokens[3]);
0070       } else {
0071         throw std::runtime_error("Invalid field map format");
0072       }
0073     }
0074     rPos.push_back(r);
0075     zPos.push_back(z);
0076     bField.push_back(Vector2(br, bz));
0077   }
0078   map_file.close();
0079   rPos.shrink_to_fit();
0080   zPos.shrink_to_fit();
0081   bField.shrink_to_fit();
0082
0083   /// [2] use helper function in core
0084   return fieldMapRZ(localToGlobalBin, rPos, zPos, bField, lengthUnit,
0085                     BFieldUnit, firstQuadrant);
0086 }
0087
0088 Acts::InterpolatedBFieldMap<
0089     Acts::Grid<Acts::Vector3, Acts::Axis<Acts::AxisType::Equidistant>,
0090                Acts::Axis<Acts::AxisType::Equidistant>,
0091                Acts::Axis<Acts::AxisType::Equidistant>>>
0092 Acts::makeMagneticFieldMapXyzFromText(
0093     const std::function<std::size_t(std::array<std::size_t, 3> binsXYZ,
0094                                     std::array<std::size_t, 3> nBinsXYZ)>&
0095         localToGlobalBin,
0096     const std::string& fieldMapFile, double lengthUnit, double BFieldUnit,
0097     bool firstOctant, const std::string& delimiter) {
0098   /// [1] Read in field map file
0099   // Grid position points in x, y and z
0100   std::vector<double> xPos;
0101   std::vector<double> yPos;
0102   std::vector<double> zPos;
0103   // components of magnetic field on grid points
0104   std::vector<Vector3> bField;
0105   // reserve estimated size
0106   xPos.reserve(kDefaultSize);
0107   yPos.reserve(kDefaultSize);
0108   zPos.reserve(kDefaultSize);
0109   bField.reserve(kDefaultSize);
0110   // [1] Read in file and fill values
0111   std::ifstream map_file(fieldMapFile.c_str(), std::ios::in);
0112   std::string line;
0113   double x = 0., y = 0., z = 0.;
0114   double bx = 0., by = 0., bz = 0.;
0115   while (std::getline(map_file, line)) {
0116     if (ignoreLine(line)) {
0117       continue;
0118     }
0119
0120     std::istringstream tmp(line);
0121     if (delimiter.empty()) {
0122       tmp >> x >> y >> z >> bx >> by >> bz;
0123     } else {
0124       std::vector<std::string> tokens;
0125       boost::split(tokens, line, boost::is_any_of(delimiter));
0126       if (tokens.size() == 6) {
0127         x = std::stod(tokens[0]);
0128         y = std::stod(tokens[1]);
0129         z = std::stod(tokens[2]);
0130         bx = std::stod(tokens[3]);
0131         by = std::stod(tokens[4]);
0132         bz = std::stod(tokens[5]);
0133       } else {
0134         throw std::runtime_error("Invalid field map format");
0135       }
0136     }
0137     xPos.push_back(x);
0138     yPos.push_back(y);
0139     zPos.push_back(z);
0140     bField.push_back(Vector3(bx, by, bz));
0141   }
0142   map_file.close();
0143   xPos.shrink_to_fit();
0144   yPos.shrink_to_fit();
0145   zPos.shrink_to_fit();
0146   bField.shrink_to_fit();
0147
0148   return fieldMapXYZ(localToGlobalBin, xPos, yPos, zPos, bField, lengthUnit,
0149                      BFieldUnit, firstOctant);
0150 }