MagFieldElement further follow ATLAS style on inline conventions

MagneticField/MagFieldElements/CMakeLists.txt

@@ -17,5 +17,4 @@ atlas_add_library( MagFieldElements
 		 
 atlas_add_test( BFieldExample_test
                 SOURCES  test/BFieldExample_test.cxx
-                INCLUDE_DIRS 
                 LINK_LIBRARIES MagFieldElements)		

MagneticField/MagFieldElements/MagFieldElements/BFieldCache.h

@@ -19,7 +19,6 @@
 
 #include "CxxUtils/restrict.h"
 #include "MagFieldElements/BFieldVector.h"
-#include <cmath>
 
 class BFieldCache
 {
@@ -28,6 +27,7 @@ public:
   BFieldCache() = default;
   // make this cache invalid, so that inside() will fail
   void invalidate();
+
   // set the z, r, phi range that defines the bin
   void setRange(double zmin,
                 double zmax,
@@ -53,9 +53,6 @@ public:
             double* ATH_RESTRICT B,
             double* ATH_RESTRICT deriv = nullptr) const;
 
-  // set the field values at each corner (rescale for current scale factor)
-  void printField() const;
-
 private:
   // bin range in z
   double m_zmin = 0.0;

MagneticField/MagFieldElements/MagFieldElements/BFieldCache.icc

 /*
   Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
-
+#include <cmath>
 inline void
 BFieldCache::invalidate()
 {
@@ -56,102 +56,4 @@ BFieldCache::inside(double z, double r, double phi) const
           r >= m_rmin && r <= m_rmax);
 }
 
-inline void
-BFieldCache::getB(const double* ATH_RESTRICT xyz,
-                  double r,
-                  double phi,
-                  double* ATH_RESTRICT B,
-                  double* ATH_RESTRICT deriv) const
-{
-
-  const double x = xyz[0];
-  const double y = xyz[1];
-  const double z = xyz[2];
-
-  // make sure phi is inside [m_phimin,m_phimax]
-  if (phi < m_phimin) {
-    phi += 2 * M_PI;
-  }
-  // fractional position inside this bin
-  const double fz = (z - m_zmin) * m_invz;
-  const double gz = 1.0 - fz;
-  const double fr = (r - m_rmin) * m_invr;
-  const double gr = 1.0 - fr;
-  const double fphi = (phi - m_phimin) * m_invphi;
-  const double gphi = 1.0 - fphi;
-  const double scale = m_scale;
-  // interpolate field values in z, r, phi
-  double Bzrphi[3];
-  for (int i = 0; i < 3; ++i) { // z, r, phi components
-    const double* field = m_field[i];
-    Bzrphi[i] = scale * (gz * (gr * (gphi * field[0] + fphi * field[1]) +
-                               fr * (gphi * field[2] + fphi * field[3])) +
-                         fz * (gr * (gphi * field[4] + fphi * field[5]) +
-                               fr * (gphi * field[6] + fphi * field[7])));
-  }
-  // convert (Bz,Br,Bphi) to (Bx,By,Bz)
-  double invr;
-  double c;
-  double s;
-  if (r > 0.0) {
-    invr = 1.0 / r;
-    c = x * invr;
-    s = y * invr;
-  } else {
-    invr = 0.0;
-    c = cos(m_phimin);
-    s = sin(m_phimin);
-  }
-  B[0] = Bzrphi[1] * c - Bzrphi[2] * s;
-  B[1] = Bzrphi[1] * s + Bzrphi[2] * c;
-  B[2] = Bzrphi[0];
-
-  // compute field derivatives if requested
-  if (deriv) {
-    const double sz = m_scale * m_invz;
-    const double sr = m_scale * m_invr;
-    const double sphi = m_scale * m_invphi;
-
-    std::array<double, 3> dBdz;
-    std::array<double, 3> dBdr;
-    std::array<double, 3> dBdphi;
-
-    for (int j = 0; j < 3; ++j) { // Bz, Br, Bphi components
-      const double* field = m_field[j];
-      dBdz[j] =
-        sz *
-        (gr * (gphi * (field[4] - field[0]) + fphi * (field[5] - field[1])) +
-         fr * (gphi * (field[6] - field[2]) + fphi * (field[7] - field[3])));
-      dBdr[j] =
-        sr *
-        (gz * (gphi * (field[2] - field[0]) + fphi * (field[3] - field[1])) +
-         fz * (gphi * (field[6] - field[4]) + fphi * (field[7] - field[5])));
-      dBdphi[j] =
-        sphi * (gz * (gr * (field[1] - field[0]) + fr * (field[3] - field[2])) +
-                fz * (gr * (field[5] - field[4]) + fr * (field[7] - field[6])));
-    }
-    // convert to cartesian coordinates
-    const double cc = c * c;
-    const double cs = c * s;
-    const double ss = s * s;
-    const double ccinvr = cc * invr;
-    const double csinvr = cs * invr;
-    const double ssinvr = ss * invr;
-    const double sinvr = s * invr;
-    const double cinvr = c * invr;
-    deriv[0] = cc * dBdr[1] - cs * dBdr[2] - csinvr * dBdphi[1] +
-               ssinvr * dBdphi[2] + sinvr * B[1];
-    deriv[1] = cs * dBdr[1] - ss * dBdr[2] + ccinvr * dBdphi[1] -
-               csinvr * dBdphi[2] - cinvr * B[1];
-    deriv[2] = c * dBdz[1] - s * dBdz[2];
-    deriv[3] = cs * dBdr[1] + cc * dBdr[2] - ssinvr * dBdphi[1] -
-               csinvr * dBdphi[2] - sinvr * B[0];
-    deriv[4] = ss * dBdr[1] + cs * dBdr[2] + csinvr * dBdphi[1] +
-               ccinvr * dBdphi[2] + cinvr * B[0];
-    deriv[5] = s * dBdz[1] + c * dBdz[2];
-    deriv[6] = c * dBdr[0] - sinvr * dBdphi[0];
-    deriv[7] = s * dBdr[0] + cinvr * dBdphi[0];
-    deriv[8] = dBdz[0];
-  }
-}
 

MagneticField/MagFieldElements/MagFieldElements/BFieldCacheZR.h

@@ -22,40 +22,16 @@ class BFieldCacheZR
 {
 public:
   // default constructor sets unphysical boundaries, so that inside() will fail
-  BFieldCacheZR()
-    : m_zmin(0.0)
-    , m_zmax(-1.0)
-    , m_rmin(0.0)
-    , m_rmax(-1.0)
-  {}
+  BFieldCacheZR();
   // invalidate this cache, so that inside() will fail
-  void invalidate()
-  {
-    m_rmin = 0.0;
-    m_rmax = -1.0;
-  }
+  void invalidate();
   // set the z, r range that defines the bin
-  void setRange(double zmin, double zmax, double rmin, double rmax)
-  {
-    m_zmin = zmin;
-    m_zmax = zmax;
-    m_rmin = rmin;
-    m_rmax = rmax;
-    m_invz = 1.0 / (zmax - zmin);
-    m_invr = 1.0 / (rmax - rmin);
-  }
+  void setRange(double zmin, double zmax, double rmin, double rmax);
   // set the field values at each corner (rescale for current scale factor)
-  void setField(int i, const BFieldVectorZR& field, double scaleFactor = 1.0)
-  {
-    m_field[0][i] = scaleFactor * field[0];
-    m_field[1][i] = scaleFactor * field[1];
-  }
+  void setField(int i, const BFieldVectorZR& field, double scaleFactor = 1.0);
   // set the multiplicative factor for the field vectors
   // test if (z, r) is inside this bin
-  bool inside(double z, double r) const
-  {
-    return (z >= m_zmin && z <= m_zmax && r >= m_rmin && r <= m_rmax);
-  }
+  bool inside(double z, double r) const;
   // interpolate the field and return B[3].
   // also compute field derivatives if deriv[9] is given.
   void getB(const double*  ATH_RESTRICT xyz,

MagneticField/MagFieldElements/MagFieldElements/BFieldCacheZR.icc

 /*
   Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
-
+inline BFieldCacheZR::BFieldCacheZR()
+  : m_zmin(0.0)
+  , m_zmax(-1.0)
+  , m_rmin(0.0)
+  , m_rmax(-1.0)
+{}
 inline void
-BFieldCacheZR::getB(const double* ATH_RESTRICT xyz,
-                    double r,
-                    double* ATH_RESTRICT B,
-                    double* ATH_RESTRICT deriv) const
+BFieldCacheZR::invalidate()
 {
-  const double x = xyz[0];
-  const double y = xyz[1];
-  const double z = xyz[2];
-  // fractional position inside this bin
-  const double fz = (z - m_zmin) * m_invz;
-  const double gz = 1.0 - fz;
-  const double fr = (r - m_rmin) * m_invr;
-  const double gr = 1.0 - fr;
-  // interpolate field values in z, r
-  double Bzr[2];
-  for (int i = 0; i < 2; ++i) { // z, r components
-    const double* field = m_field[i];
-    Bzr[i] = gz * (gr * field[0] + fr * field[1]) +
-             fz * (gr * field[2] + fr * field[3]);
-  }
-  // convert (Bz,Br) to (Bx,By,Bz)
-  double invr;
-  if (r > 0.0) {
-    invr = 1.0 / r;
-  } else {
-    invr = 0.0;
-  }
-  const double c(x * invr);
-  const double s(y * invr);
-  B[0] = Bzr[1] * c;
-  B[1] = Bzr[1] * s;
-  B[2] = Bzr[0];
-
-  // compute field derivatives if requested
-  if (deriv) {
-    std::array<double, 2> dBdz;
-    std::array<double, 2> dBdr;
-    for (int j = 0; j < 2; ++j) { // Bz, Br components
-      const double* field = m_field[j];
-      dBdz[j] =
-          m_invz * (gr * (field[2] - field[0]) +
-                    fr * (field[3] - field[1]));
-      dBdr[j] =
-          m_invr * (gz * (field[1] - field[0]) +
-                    fz * (field[3] - field[2]));
-    }
-    // convert to cartesian coordinates
-    const double cc = c * c;
-    const double cs = c * s;
-    const double ss = s * s;
-    const double sinvr = s * invr;
-    const double cinvr = c * invr;
-    deriv[0] = cc * dBdr[1] + sinvr * B[1];
-    deriv[1] = cs * dBdr[1] - cinvr * B[1];
-    deriv[2] = c * dBdz[1];
-    deriv[3] = cs * dBdr[1] - sinvr * B[0];
-    deriv[4] = ss * dBdr[1] + cinvr * B[0];
-    deriv[5] = s * dBdz[1];
-    deriv[6] = c * dBdr[0];
-    deriv[7] = s * dBdr[0];
-    deriv[8] = dBdz[0];
-  }
+  m_rmin = 0.0;
+  m_rmax = -1.0;
+}
+inline void
+BFieldCacheZR::setRange(double zmin, double zmax, double rmin, double rmax)
+{
+  m_zmin = zmin;
+  m_zmax = zmax;
+  m_rmin = rmin;
+  m_rmax = rmax;
+  m_invz = 1.0 / (zmax - zmin);
+  m_invr = 1.0 / (rmax - rmin);
+}
+// set the field values at each corner (rescale for current scale factor)
+inline void
+BFieldCacheZR::setField(int i,
+                        const BFieldVectorZR& field,
+                        double scaleFactor)
+{
+  m_field[0][i] = scaleFactor * field[0];
+  m_field[1][i] = scaleFactor * field[1];
+}
+inline bool
+BFieldCacheZR::inside(double z, double r) const
+{
+  return (z >= m_zmin && z <= m_zmax && r >= m_rmin && r <= m_rmax);
 }
 
 

MagneticField/MagFieldElements/src/BFieldCache.cxx

@@ -3,21 +3,105 @@
 */
 
 #include "MagFieldElements/BFieldCache.h"
-#include <iostream>
+#include <cmath>
 
-// set the field values at each corner (rescale for current scale factor)
 void
-BFieldCache::printField() const
+BFieldCache::getB(const double* ATH_RESTRICT xyz,
+                  double r,
+                  double phi,
+                  double* ATH_RESTRICT B,
+                  double* ATH_RESTRICT deriv) const
 {
-  // print out field values
-  std::cout << "Field at corner i, for each component j (Bz, Br, Bphi)"
-            << '\n';
 
-  for (int i = 0; i < 8; ++i) {
-    for (int j = 0; j < 3; ++j) {
-      std::cout << i << "," << j << ": " << m_field[j][i] << ", ";
+  const double x = xyz[0];
+  const double y = xyz[1];
+  const double z = xyz[2];
+
+  // make sure phi is inside [m_phimin,m_phimax]
+  if (phi < m_phimin) {
+    phi += 2 * M_PI;
+  }
+  // fractional position inside this bin
+  const double fz = (z - m_zmin) * m_invz;
+  const double gz = 1.0 - fz;
+  const double fr = (r - m_rmin) * m_invr;
+  const double gr = 1.0 - fr;
+  const double fphi = (phi - m_phimin) * m_invphi;
+  const double gphi = 1.0 - fphi;
+  const double scale = m_scale;
+  // interpolate field values in z, r, phi
+  double Bzrphi[3];
+  for (int i = 0; i < 3; ++i) { // z, r, phi components
+    const double* field = m_field[i];
+    Bzrphi[i] = scale * (gz * (gr * (gphi * field[0] + fphi * field[1]) +
+                               fr * (gphi * field[2] + fphi * field[3])) +
+                         fz * (gr * (gphi * field[4] + fphi * field[5]) +
+                               fr * (gphi * field[6] + fphi * field[7])));
+  }
+  // convert (Bz,Br,Bphi) to (Bx,By,Bz)
+  double invr;
+  double c;
+  double s;
+  if (r > 0.0) {
+    invr = 1.0 / r;
+    c = x * invr;
+    s = y * invr;
+  } else {
+    invr = 0.0;
+    c = cos(m_phimin);
+    s = sin(m_phimin);
+  }
+  B[0] = Bzrphi[1] * c - Bzrphi[2] * s;
+  B[1] = Bzrphi[1] * s + Bzrphi[2] * c;
+  B[2] = Bzrphi[0];
+
+  // compute field derivatives if requested
+  if (deriv) {
+    const double sz = m_scale * m_invz;
+    const double sr = m_scale * m_invr;
+    const double sphi = m_scale * m_invphi;
+
+    std::array<double, 3> dBdz;
+    std::array<double, 3> dBdr;
+    std::array<double, 3> dBdphi;
+
+    for (int j = 0; j < 3; ++j) { // Bz, Br, Bphi components
+      const double* field = m_field[j];
+      dBdz[j] =
+        sz *
+        (gr * (gphi * (field[4] - field[0]) + fphi * (field[5] - field[1])) +
+         fr * (gphi * (field[6] - field[2]) + fphi * (field[7] - field[3])));
+      dBdr[j] =
+        sr *
+        (gz * (gphi * (field[2] - field[0]) + fphi * (field[3] - field[1])) +
+         fz * (gphi * (field[6] - field[4]) + fphi * (field[7] - field[5])));
+      dBdphi[j] =
+        sphi * (gz * (gr * (field[1] - field[0]) + fr * (field[3] - field[2])) +
+                fz * (gr * (field[5] - field[4]) + fr * (field[7] - field[6])));
     }
-    std::cout << '\n';
+    // convert to cartesian coordinates
+    const double cc = c * c;
+    const double cs = c * s;
+    const double ss = s * s;
+    const double ccinvr = cc * invr;
+    const double csinvr = cs * invr;
+    const double ssinvr = ss * invr;
+    const double sinvr = s * invr;
+    const double cinvr = c * invr;
+    deriv[0] = cc * dBdr[1] - cs * dBdr[2] - csinvr * dBdphi[1] +
+               ssinvr * dBdphi[2] + sinvr * B[1];
+    deriv[1] = cs * dBdr[1] - ss * dBdr[2] + ccinvr * dBdphi[1] -
+               csinvr * dBdphi[2] - cinvr * B[1];
+    deriv[2] = c * dBdz[1] - s * dBdz[2];
+    deriv[3] = cs * dBdr[1] + cc * dBdr[2] - ssinvr * dBdphi[1] -
+               csinvr * dBdphi[2] - sinvr * B[0];
+    deriv[4] = ss * dBdr[1] + cs * dBdr[2] + csinvr * dBdphi[1] +
+               ccinvr * dBdphi[2] + cinvr * B[0];
+    deriv[5] = s * dBdz[1] + c * dBdz[2];
+    deriv[6] = c * dBdr[0] - sinvr * dBdphi[0];
+    deriv[7] = s * dBdr[0] + cinvr * dBdphi[0];
+    deriv[8] = dBdz[0];
   }
 }
 
+

MagneticField/MagFieldElements/src/BFieldCacheZR.cxx

+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "MagFieldElements/BFieldCacheZR.h"
+#include <cmath>
+
+void
+BFieldCacheZR::getB(const double* ATH_RESTRICT xyz,
+                    double r,
+                    double* ATH_RESTRICT B,
+                    double* ATH_RESTRICT deriv) const
+{
+  const double x = xyz[0];
+  const double y = xyz[1];
+  const double z = xyz[2];
+  // fractional position inside this bin
+  const double fz = (z - m_zmin) * m_invz;
+  const double gz = 1.0 - fz;
+  const double fr = (r - m_rmin) * m_invr;
+  const double gr = 1.0 - fr;
+  // interpolate field values in z, r
+  double Bzr[2];
+  for (int i = 0; i < 2; ++i) { // z, r components
+    const double* field = m_field[i];
+    Bzr[i] = gz * (gr * field[0] + fr * field[1]) +
+             fz * (gr * field[2] + fr * field[3]);
+  }
+  // convert (Bz,Br) to (Bx,By,Bz)
+  double invr;
+  if (r > 0.0) {
+    invr = 1.0 / r;
+  } else {
+    invr = 0.0;
+  }
+  const double c(x * invr);
+  const double s(y * invr);
+  B[0] = Bzr[1] * c;
+  B[1] = Bzr[1] * s;
+  B[2] = Bzr[0];
+
+  // compute field derivatives if requested
+  if (deriv) {
+    std::array<double, 2> dBdz;
+    std::array<double, 2> dBdr;
+    for (int j = 0; j < 2; ++j) { // Bz, Br components
+      const double* field = m_field[j];
+      dBdz[j] =
+        m_invz * (gr * (field[2] - field[0]) + fr * (field[3] - field[1]));
+      dBdr[j] =
+        m_invr * (gz * (field[1] - field[0]) + fz * (field[3] - field[2]));
+    }
+    // convert to cartesian coordinates
+    const double cc = c * c;
+    const double cs = c * s;
+    const double ss = s * s;
+    const double sinvr = s * invr;
+    const double cinvr = c * invr;
+    deriv[0] = cc * dBdr[1] + sinvr * B[1];
+    deriv[1] = cs * dBdr[1] - cinvr * B[1];
+    deriv[2] = c * dBdz[1];
+    deriv[3] = cs * dBdr[1] - sinvr * B[0];
+    deriv[4] = ss * dBdr[1] + cinvr * B[0];
+    deriv[5] = s * dBdz[1];
+    deriv[6] = c * dBdr[0];
+    deriv[7] = s * dBdr[0];
+    deriv[8] = dBdz[0];
+  }
+}
+

MagneticField/MagFieldElements/test/BFieldExample_test.cxx

@@ -122,15 +122,6 @@ public:
     BFieldCache cache3d;
     double z{ 0 }, r{ 1250 }, phi{ 1.6 };
 
-    // if (doDebug) std::cout << "do getCache old " << '\n';
-    // zone.m_doNew = false;
-    // zone.getCache(z, r, phi, cache3d, 1);
-    // cache3d.printField();
-    // if (doDebug) std::cout << "do getCache new " << '\n';
-    // zone.m_doNew = true;
-    // zone.getCache(z, r, phi, cache3d, 1);
-    // cache3d.printField();
-
     // get field at steps of 10 mm from 1200 to 1300
     int status{ 0 };