From 930da18c49e5dc68b6d6db0859209edbf67200d7 Mon Sep 17 00:00:00 2001 From: christos <christos@cern.ch> Date: Fri, 4 Sep 2020 03:46:35 +0100 Subject: [PATCH] MagFieldElement move small inline methods to .icc files, move a bit larger methods to .cxx --- .../MagFieldElements/AtlasFieldCache.h | 28 +++-- .../MagFieldElements/AtlasFieldCache.icc | 116 +++--------------- .../MagFieldElements/BFieldZone.h | 49 ++------ .../MagFieldElements/BFieldZone.icc | 49 ++++++++ .../MagFieldElements/src/AtlasFieldCache.cxx | 109 +++++++++++++--- 5 files changed, 186 insertions(+), 165 deletions(-) create mode 100644 MagneticField/MagFieldElements/MagFieldElements/BFieldZone.icc diff --git a/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.h b/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.h index 169a0140221..7161f00aa43 100644 --- a/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.h +++ b/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.h @@ -37,9 +37,10 @@ namespace MagField { class AtlasFieldCache { public: - AtlasFieldCache(); - // ** constructor to setup with field scale and magnetic field service for - // first access to field */ + AtlasFieldCache() = default; + /** constructor to setup with field scale + * and magnetic field service for + * first access to field */ AtlasFieldCache(double solFieldScale, double torFieldScale, const AtlasFieldMap* fieldMap); @@ -54,16 +55,17 @@ public: /// Temporary flag for switching between 'old' and 'new' magField usage bool useNewBfieldCache() { return true; } - /** get B field value at given position */ - /** xyz[3] is in mm, bxyz[3] is in kT */ - /** if deriv[9] is given, field derivatives are returned in kT/mm */ - inline void getField(const double* ATH_RESTRICT xyz, - double* ATH_RESTRICT bxyz, - double* ATH_RESTRICT deriv = nullptr); + /** get B field value at given position + * xyz[3] is in mm, bxyz[3] is in kT + * if deriv[9] is given, field derivatives are returned in kT/mm + * */ + void getField(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT bxyz, + double* ATH_RESTRICT deriv = nullptr); - inline void getFieldZR(const double* ATH_RESTRICT xyz, - double* ATH_RESTRICT bxyz, - double* ATH_RESTRICT deriv = nullptr); + void getFieldZR(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT bxyz, + double* ATH_RESTRICT deriv = nullptr); /** status of the magnets */ bool solenoidOn() const; @@ -73,7 +75,9 @@ private: AtlasFieldCache(const AtlasFieldCache& other) = delete; AtlasFieldCache& operator=(const AtlasFieldCache& other) = delete; + /// fill given magnetic field zone */ bool fillFieldCache(double z, double r, double phi); + /// fill Z-R cache for solenoid */ bool fillFieldCacheZR(double z, double r); /// Full 3d field diff --git a/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.icc b/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.icc index 8f5bd59d021..cf32e19190d 100644 --- a/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.icc +++ b/MagneticField/MagFieldElements/MagFieldElements/AtlasFieldCache.icc @@ -1,8 +1,24 @@ /* Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration */ +inline MagField::AtlasFieldCache::AtlasFieldCache(double solFieldScale, + double torFieldScale, + const AtlasFieldMap* fieldMap) + : m_solScale(solFieldScale) + , m_torScale(torFieldScale) + , + // set field service + m_fieldMap(fieldMap) + +{ + if (m_fieldMap) { + // save ZR bfield + m_meshZR = m_fieldMap->getBFieldMesh(); + // Get solenoid zone id from field service + m_solZoneId = fieldMap->solenoidZoneId(); + } +} -/** fill given magnetic field zone */ inline bool MagField::AtlasFieldCache::fillFieldCache(double z, double r, double phi) { @@ -27,7 +43,6 @@ MagField::AtlasFieldCache::fillFieldCache(double z, double r, double phi) return true; } -/** fill Z-R cache for solenoid */ inline bool MagField::AtlasFieldCache::fillFieldCacheZR(double z, double r) { @@ -44,103 +59,6 @@ MagField::AtlasFieldCache::fillFieldCacheZR(double z, double r) return true; } -inline void -MagField::AtlasFieldCache::getField(const double* ATH_RESTRICT xyz, - double* ATH_RESTRICT bxyz, - double* ATH_RESTRICT deriv) -{ - // Allow for the case of no map for testing - if (m_fieldMap == nullptr) { - // return 0 bfield if map is missing - bxyz[0] = bxyz[1] = bxyz[2] = 0; - if (deriv) { - for (int i = 0; i < 9; i++) { - deriv[i] = 0.; - } - } - return; - } - - const double x = xyz[0]; - const double y = xyz[1]; - const double z = xyz[2]; - const double r = std::sqrt(x * x + y * y); - const double phi = std::atan2(y, x); - - // test if initialised and the cache is valid - if (!m_cache3d.inside(z, r, phi)) { - // cache is invalid -> refresh cache - if (!fillFieldCache(z, r, phi)) { - // caching failed -> outside the valid map volume - // return default field (0.1 gauss) - const double defaultB(0.1 * Gaudi::Units::gauss); - bxyz[0] = bxyz[1] = bxyz[2] = defaultB; - // return zero gradient if requested - if (deriv) { - for (int i = 0; i < 9; i++) { - deriv[i] = 0.; - } - } - return; - } - } - - // do interpolation (cache3d has correct scale factor) - m_cache3d.getB(xyz, r, phi, bxyz, deriv); - - if (!m_cond) { - return; - } - // add biot savart component - must add in scale factor to avoid changing - // conductor SF since the conductor is part of the static magnetic field model - const size_t condSize = m_cond->size(); - for (size_t i = 0; i < condSize; i++) { - (*m_cond)[i].addBiotSavart(m_scaleToUse, xyz, bxyz, deriv); - } -} - -inline void -MagField::AtlasFieldCache::getFieldZR(const double* ATH_RESTRICT xyz, - double* ATH_RESTRICT bxyz, - double* ATH_RESTRICT deriv) -{ - - // Allow for the case of no map for testing - if (m_fieldMap == nullptr) { - // constant ATLAS magnetic field if no map has been set - for testing - constexpr double TEST_BFIELD = 1.997; - bxyz[0] = bxyz[1] = 0; - bxyz[2] = TEST_BFIELD; - if (deriv) { - for (int i = 0; i < 9; i++) { - deriv[i] = 0.; - } - } - return; - } - - const double x = xyz[0]; - const double y = xyz[1]; - const double z = xyz[2]; - const double r = sqrt(x * x + y * y); - - // test if the cache was initialized and the ZR cache is valid for current - // position - if (!m_cacheZR.inside(z, r)) { - - // cache is invalid -> refresh cache - if (!fillFieldCacheZR(z, r)) { - - // caching failed -> outside the valid z-r map volume - // call the full version of getField() - getField(xyz, bxyz, deriv); - return; - } - } - // do interpolation - m_cacheZR.getB(xyz, r, bxyz, deriv); -} - inline bool MagField::AtlasFieldCache::solenoidOn() const { diff --git a/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.h b/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.h index baa556d12ae..57fbb3ea623 100644 --- a/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.h +++ b/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.h @@ -9,9 +9,11 @@ // // Masahiro Morii, Harvard University // +// #ifndef BFIELDZONE_H #define BFIELDZONE_H +#include "CxxUtils/restrict.h" #include "MagFieldElements/BFieldCond.h" #include "MagFieldElements/BFieldMesh.h" #include <vector> @@ -34,58 +36,29 @@ public: ; } // add elements to vectors - void appendCond(const BFieldCond& cond) { m_cond.push_back(cond); } + void appendCond(const BFieldCond& cond); // compute Biot-Savart magnetic field and add to B[3] - inline void addBiotSavart(const double* xyz, - double* B, - double* deriv = nullptr) const; + void addBiotSavart(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT B, + double* ATH_RESTRICT deriv = nullptr) const; // scale B field by a multiplicative factor: RDS 2019/09 - no longer used. // Scaling is done in cachec - void scaleField(double factor) - { - scaleBscale(factor); - for (unsigned i = 0; i < ncond(); i++) { - m_cond[i].scaleCurrent(factor); - } - } + void scaleField(double factor); // accessors int id() const { return m_id; } unsigned ncond() const { return m_cond.size(); } const BFieldCond& cond(int i) const { return m_cond[i]; } const std::vector<BFieldCond>* condVector() const { return &m_cond; } - int memSize() const - { - return BFieldMesh<short>::memSize() + sizeof(int) + - sizeof(BFieldCond) * m_cond.capacity(); - } + int memSize() const; // adjust the min/max edges to a new value - void adjustMin(int i, double x) - { - m_min[i] = x; - m_mesh[i].front() = x; - } - void adjustMax(int i, double x) - { - m_max[i] = x; - m_mesh[i].back() = x; - } + void adjustMin(int i, double x); + void adjustMax(int i, double x); private: int m_id; // zone ID number std::vector<BFieldCond> m_cond; // list of current conductors }; -// inline functions - -// -// Compute magnetic field due to the conductors -// -void -BFieldZone::addBiotSavart(const double* xyz, double* B, double* deriv) const -{ - for (unsigned i = 0; i < m_cond.size(); i++) { - m_cond[i].addBiotSavart(1, xyz, B, deriv); - } -} +#include "BFieldZone.icc" #endif diff --git a/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.icc b/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.icc new file mode 100644 index 00000000000..071400e3729 --- /dev/null +++ b/MagneticField/MagFieldElements/MagFieldElements/BFieldZone.icc @@ -0,0 +1,49 @@ +/* + Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration +*/ +inline void +BFieldZone::appendCond(const BFieldCond& cond) +{ + m_cond.push_back(cond); +} + +inline void +BFieldZone::scaleField(double factor) +{ + scaleBscale(factor); + for (unsigned i = 0; i < ncond(); i++) { + m_cond[i].scaleCurrent(factor); + } +} + +inline void +BFieldZone::addBiotSavart(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT B, + double* ATH_RESTRICT deriv) const +{ + for (unsigned i = 0; i < m_cond.size(); i++) { + m_cond[i].addBiotSavart(1, xyz, B, deriv); + } +} + +inline int +BFieldZone::memSize() const +{ + return BFieldMesh<short>::memSize() + sizeof(int) + + sizeof(BFieldCond) * m_cond.capacity(); +} + +inline void +BFieldZone::adjustMin(int i, double x) +{ + m_min[i] = x; + m_mesh[i].front() = x; +} + +inline void +BFieldZone::adjustMax(int i, double x) +{ + m_max[i] = x; + m_mesh[i].back() = x; +} + diff --git a/MagneticField/MagFieldElements/src/AtlasFieldCache.cxx b/MagneticField/MagFieldElements/src/AtlasFieldCache.cxx index 89ace3a1077..a3927c4553d 100644 --- a/MagneticField/MagFieldElements/src/AtlasFieldCache.cxx +++ b/MagneticField/MagFieldElements/src/AtlasFieldCache.cxx @@ -11,23 +11,100 @@ // #include "MagFieldElements/AtlasFieldCache.h" -/// Constructor -MagField::AtlasFieldCache::AtlasFieldCache() = default; - -MagField::AtlasFieldCache::AtlasFieldCache(double solFieldScale, - double torFieldScale, - const AtlasFieldMap* fieldMap) - : - m_solScale(solFieldScale), - m_torScale(torFieldScale), - // set field service - m_fieldMap(fieldMap) +void +MagField::AtlasFieldCache::getField(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT bxyz, + double* ATH_RESTRICT deriv) +{ + // Allow for the case of no map for testing + if (m_fieldMap == nullptr) { + // return 0 bfield if map is missing + bxyz[0] = bxyz[1] = bxyz[2] = 0; + if (deriv) { + for (int i = 0; i < 9; i++) { + deriv[i] = 0.; + } + } + return; + } + + const double x = xyz[0]; + const double y = xyz[1]; + const double z = xyz[2]; + const double r = std::sqrt(x * x + y * y); + const double phi = std::atan2(y, x); + + // test if initialised and the cache is valid + if (!m_cache3d.inside(z, r, phi)) { + // cache is invalid -> refresh cache + if (!fillFieldCache(z, r, phi)) { + // caching failed -> outside the valid map volume + // return default field (0.1 gauss) + const double defaultB(0.1 * Gaudi::Units::gauss); + bxyz[0] = bxyz[1] = bxyz[2] = defaultB; + // return zero gradient if requested + if (deriv) { + for (int i = 0; i < 9; i++) { + deriv[i] = 0.; + } + } + return; + } + } + + // do interpolation (cache3d has correct scale factor) + m_cache3d.getB(xyz, r, phi, bxyz, deriv); + + if (!m_cond) { + return; + } + // add biot savart component - must add in scale factor to avoid changing + // conductor SF since the conductor is part of the static magnetic field model + const size_t condSize = m_cond->size(); + for (size_t i = 0; i < condSize; i++) { + (*m_cond)[i].addBiotSavart(m_scaleToUse, xyz, bxyz, deriv); + } +} +void +MagField::AtlasFieldCache::getFieldZR(const double* ATH_RESTRICT xyz, + double* ATH_RESTRICT bxyz, + double* ATH_RESTRICT deriv) { - if (m_fieldMap) { - // save ZR bfield - m_meshZR = m_fieldMap->getBFieldMesh() ; - // Get solenoid zone id from field service - m_solZoneId = fieldMap->solenoidZoneId(); + + // Allow for the case of no map for testing + if (m_fieldMap == nullptr) { + // constant ATLAS magnetic field if no map has been set - for testing + constexpr double TEST_BFIELD = 1.997; + bxyz[0] = bxyz[1] = 0; + bxyz[2] = TEST_BFIELD; + if (deriv) { + for (int i = 0; i < 9; i++) { + deriv[i] = 0.; + } } + return; + } + + const double x = xyz[0]; + const double y = xyz[1]; + const double z = xyz[2]; + const double r = sqrt(x * x + y * y); + + // test if the cache was initialized and the ZR cache is valid for current + // position + if (!m_cacheZR.inside(z, r)) { + + // cache is invalid -> refresh cache + if (!fillFieldCacheZR(z, r)) { + + // caching failed -> outside the valid z-r map volume + // call the full version of getField() + getField(xyz, bxyz, deriv); + return; + } + } + // do interpolation + m_cacheZR.getB(xyz, r, bxyz, deriv); } + -- GitLab