diff --git a/Detector/MagneticFieldMap/CMakeLists.txt b/Detector/MagneticFieldMap/CMakeLists.txt index ac154c8ae17f339802d26b7a1e5f638a87c62542..7c65b50318a2d0040922057a2fcc36a84db3ec50 100644 --- a/Detector/MagneticFieldMap/CMakeLists.txt +++ b/Detector/MagneticFieldMap/CMakeLists.txt @@ -5,6 +5,7 @@ gaudi_add_module(MagneticFieldMap SOURCES src/GenericBFieldMapBrBz.cpp src/GenericBFieldMapBrBzFactory.cpp + src/FieldMapFileProvider.cpp LINK Gaudi::GaudiKernel ${DD4hep_COMPONENT_LIBRARIES} ${ROOT_LIBRARIES} diff --git a/Detector/MagneticFieldMap/src/FieldMapFileProvider.cpp b/Detector/MagneticFieldMap/src/FieldMapFileProvider.cpp new file mode 100644 index 0000000000000000000000000000000000000000..388ea22e00e2ddeb218178a7330a92702472413e --- /dev/null +++ b/Detector/MagneticFieldMap/src/FieldMapFileProvider.cpp @@ -0,0 +1,229 @@ +#include "FieldMapFileProvider.h" + +#include <cmath> +#include <string> +#include <map> +#include <vector> +#include <iostream> +#include <fstream> +#include <sstream> + +FieldMapFileProvider::FieldMapFileProvider(const std::string& url_) + : m_url(url_), rBinMax(-1), zBinMax(-1), drBin(-1), dzBin(-1) { + init(); +} + +int FieldMapFileProvider::rBinIdx(double r) { + double absr = std::fabs(r); + + return -1; +} + +int FieldMapFileProvider::zBinIdx(double z) { + double absz = std::fabs(z); + + return -1; +} + +void FieldMapFileProvider::access(int rbin, int zbin, double& Br, double& Bz) { + +} + +// ====================== +// Below are private impl +// ====================== + +void FieldMapFileProvider::init() { + // parse the url + // Br=/tmp/lint/CEPCSW/Br.csv;Bz=/tmp/lint/CEPCSW/Bz.csv + + std::map<std::string, std::string> url_parsed; + + size_t idx_begin = 0; + size_t idx_end = m_url.find(";"); + + while (true) { + std::string keyval = m_url.substr(idx_begin, idx_end-idx_begin); + + if (keyval.size()) { + std::cout << "---> keyval: " << keyval << std::endl; + // separate by '=' + size_t idx_sep = keyval.find("="); + if (idx_sep == std::string::npos) { + std::string error_msg = "[ERROR] FieldMapFileProvider: Pleaes specify label=path. "; + throw std::runtime_error(error_msg); + + } + std::string key = keyval.substr(0, idx_sep); + std::string val = keyval.substr(idx_sep+1); + std::cout << "-----> key: " << key << std::endl; + std::cout << "-----> val: " << val << std::endl; + // insert into map + + // duplicated + if (url_parsed.count(key)) { + std::string error_msg = "[ERROR] FieldMapFileProvider: duplicated key '" + key + "'"; + throw std::runtime_error(error_msg); + } + + url_parsed[key] = val; + } + + if (idx_end == std::string::npos) { + break; + } + + idx_begin = idx_end+1; + idx_end = m_url.find(';', idx_begin); + } + + // access the map 'url_parsed' + std::string key; + + // load Br + key = "Br"; + if (url_parsed.count(key) == 0) { + std::string error_msg = "[ERROR] FieldMapFileProvider: missing key '" + key + "'"; + throw std::runtime_error(error_msg); + } + loadBr(url_parsed[key]); + + // load Bz + key = "Bz"; + if (url_parsed.count(key) == 0) { + std::string error_msg = "[ERROR] FieldMapFileProvider: missing key '" + key + "'"; + throw std::runtime_error(error_msg); + } + loadBz(url_parsed[key]); + +} + + + +bool FieldMapFileProvider::loadBr(const std::string& fn) { + + // ----> r + // | + // v + // z + // the shape + int ncol = -1; + int nrow = -1; + + // load data + bool st = loadCSV(fn, Brdata, ncol, nrow); + + // update the metadata + // for (int i = 0; i < ncol; ++i) { + // std::cout << "i: " << i << " data[i]: " << data[i] << std::endl; + // } + nr = ncol-1; // skip the top row and left col + nz = nrow-1; + + rBinMax = Brdata[ncol-1]; + zBinMax = Brdata[ncol*(nrow-1)]; + + drBin = rBinMax / (nr-1); + dzBin = zBinMax / (nz-1); + + std::cout << "nr: " << nr << std::endl; + std::cout << "nz: " << nz << std::endl; + std::cout << "rBinMax: " << rBinMax << std::endl; + std::cout << "zBinMax: " << zBinMax << std::endl; + std::cout << "drBin: " << drBin << std::endl; + std::cout << "dzBin: " << dzBin << std::endl; + + + + return true; +} + +bool FieldMapFileProvider::loadBz(const std::string& fn) { + + // the shape + int ncol = -1; + int nrow = -1; + bool st = loadCSV(fn, Bzdata, ncol, nrow); + + return true; +} + +bool FieldMapFileProvider::loadCSV(const std::string& fn, + std::vector<double>& data, + int& ncol, int& nrow) { + + std::ifstream input(fn); + std::string tmpline; + + ncol = 0; + nrow = 0; + + + while(std::getline(input, tmpline)) { + std::cout << "------> " << tmpline << std::endl; + ++nrow; + + // split by , + int col_per_line = 0; + + size_t idx_begin = 0; + size_t idx_end = tmpline.find(","); + + while (true) { + std::cout << "----------------> idx: " + << idx_begin << "->" << idx_end << std::endl; + std::string val = tmpline.substr(idx_begin, idx_end-idx_begin); + double v = 0.0; + // parse val + std::stringstream ss; + ss << val; + ss >> v; + // if (not ss.good()) { + // v = 0; + // ss.clear(); + // } + + std::cout << "----------------> v: " << v << std::endl; + // if (val.size()==0 || ) { // if it is empty or not number, use 0 + // data.push_back(0.0); + // } else { + // } + data.push_back(v); + + ++col_per_line; + + // break + if (idx_end == std::string::npos) { + break; + } + + // goto next val + idx_begin = idx_end+1; + idx_end = tmpline.find(',', idx_begin); + + } + + // if it is the first line, we need to store the column sizes + if (ncol == 0) { + ncol = col_per_line; + } + + std::cout << "--> ncol: " << ncol + << " col_per_line: " << col_per_line + << std::endl; + + if (ncol != col_per_line) { + std::string error_msg = "[ERROR] FieldMapFileProvider: Mismatch columns. "; + throw std::runtime_error(error_msg); + } + + } + + std::cout << "SUMMARY: " + << " size of the data: " << data.size() + << " shape (ncol: " << ncol + << " , nrow: " << nrow + << ")" << std::endl; + + return true; +} diff --git a/Detector/MagneticFieldMap/src/FieldMapFileProvider.h b/Detector/MagneticFieldMap/src/FieldMapFileProvider.h new file mode 100644 index 0000000000000000000000000000000000000000..3661863bf49032ebe7c6b1f0c5943c153b856ad7 --- /dev/null +++ b/Detector/MagneticFieldMap/src/FieldMapFileProvider.h @@ -0,0 +1,49 @@ +#ifndef FieldMapFileProvider_h +#define FieldMapFileProvider_h + +#include "IFieldMapProvider.h" + +#include <string> +#include <vector> + +class FieldMapFileProvider: public IFieldMapProvider { + +public: + FieldMapFileProvider(const std::string& url_); + + // Meta data about the map + virtual int rBinIdx(double r); + virtual int zBinIdx(double z); + + // The Br and Bz + virtual void access(int rbin, int zbin, double& Br, double& Bz); + +private: + void init(); + + bool loadBr(const std::string& fn); + bool loadBz(const std::string& fn); + + bool loadCSV(const std::string& fn, std::vector<double>& data, int& ncol, int& nrow); +private: + std::string m_url; // the url passed from GenericBFieldMapBrBz + + std::string m_filename_Br; + std::string m_filename_Bz; + + // the data will include the left col (z) and top row (r) + std::vector<double> Brdata; + std::vector<double> Bzdata; + + // in this impl, the rBin/zBin should be consistent between Br and Bz. + int nr; // include the two endpoints, so nr = nrBin + 1 + int nz; + + double rBinMax; + double zBinMax; + + double drBin; + double dzBin; +}; + +#endif diff --git a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.cpp b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.cpp index fe268217c0b0376cd09c45ff3d3cf51db36f7585..568d70c80d14275e7bd2bd6d41b868e9844c147a 100644 --- a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.cpp +++ b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.cpp @@ -1,6 +1,8 @@ #include "GenericBFieldMapBrBz.h" +#include "FieldMapFileProvider.h" + GenericBFieldMapBrBz::GenericBFieldMapBrBz() : m_provider(nullptr) { type = dd4hep::CartesianField::MAGNETIC; @@ -17,9 +19,10 @@ void GenericBFieldMapBrBz::fieldComponents(const double* pos, double* field) { return; } -void GenericBFieldMapBrBz::init_provider(const std::string& provider) { +void GenericBFieldMapBrBz::init_provider(const std::string& provider, const std::string& url) { if (provider == "file") { std::cout << "Initialize provider with file. " << std::endl; + m_provider = new FieldMapFileProvider(url); } else if (provider == "db") { std::cout << "Initialize provider with file. " << std::endl; } else { diff --git a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.h b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.h index 3b1711f5abcf411253dfd8901510510fab8bfd0c..93abfcfbbb5650d5e5945b774280e6586b1c17c3 100644 --- a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.h +++ b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBz.h @@ -24,7 +24,7 @@ public: public: // following are interfaces to configure this field map - void init_provider(const std::string& provider); + void init_provider(const std::string& provider, const std::string& url); private: diff --git a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBzFactory.cpp b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBzFactory.cpp index b738332c12f8f8396aa79de1cbf9751ae906885a..d7d9fe302d75e214e290f52ebf2382441df305c0 100644 --- a/Detector/MagneticFieldMap/src/GenericBFieldMapBrBzFactory.cpp +++ b/Detector/MagneticFieldMap/src/GenericBFieldMapBrBzFactory.cpp @@ -53,11 +53,22 @@ static dd4hep::Ref_t create_GenericBFieldMapBrBz(dd4hep::Detector& , } + dd4hep::xml::Component source(xmlParameter.child(_Unicode(source))); + + // get the url + bool hasUrl = source.hasAttr(_Unicode(url)); + if (!hasUrl) { + std::string error_msg = "[ERROR] GenericBFieldMapBrBz: Must specify the 'url' in 'source' tag. "; + throw std::runtime_error(error_msg); + } + + std::string url = source.attr<std::string>(_Unicode(url)); + // 2. create the CartesianField dd4hep::CartesianField obj; GenericBFieldMapBrBz* ptr = new GenericBFieldMapBrBz(); - ptr->init_provider(provider); + ptr->init_provider(provider, url); obj.assign(ptr, xmlParameter.nameStr(), xmlParameter.typeStr()); diff --git a/Detector/MagneticFieldMap/src/IFieldMapProvider.h b/Detector/MagneticFieldMap/src/IFieldMapProvider.h index e3eda83f9653e33f8fc24967a32685b16d250d09..3d2384224d27c05c5cdd22175ef0f3f871151683 100644 --- a/Detector/MagneticFieldMap/src/IFieldMapProvider.h +++ b/Detector/MagneticFieldMap/src/IFieldMapProvider.h @@ -11,7 +11,7 @@ */ class IFieldMapProvider { - +public: // Meta data about the map virtual int rBinIdx(double r) = 0; virtual int zBinIdx(double z) = 0;