CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MFAnalyzer.cc
Go to the documentation of this file.
3 
6 
8 
12 
17 
19 {
20 public:
21  explicit MFAnalyzer( const edm::ParameterSet& );
22  ~MFAnalyzer( void );
23 
24  bool valid (void) { return m_valid; }
25 
26 private:
27  virtual void analyze( const edm::Event&, const edm::EventSetup& );
28  void evaluate( const double point[3], double field[3] ) const;
29  unsigned m_mapDensityX;
30  unsigned m_mapDensityY;
31  unsigned m_mapDensityZ;
32  double m_minX;
33  double m_maxX;
34  double m_minY;
35  double m_maxY;
36  double m_minZ;
37  double m_maxZ;
38  double m_xBaseDir;
39  double m_yBaseDir;
40  double m_zBaseDir;
41  bool m_valid;
43 };
44 
46  : m_valid( false )
47 {
48  m_mapDensityX = iPset.getUntrackedParameter<unsigned>( "mapDensityX", 10 );
49  m_mapDensityY = iPset.getUntrackedParameter<unsigned>( "mapDensityY", 10 );
50  m_mapDensityZ = iPset.getUntrackedParameter<unsigned>( "mapDensityY", 10 );
51  m_minX = iPset.getUntrackedParameter<double>( "minX", -18.0 );
52  m_maxX = iPset.getUntrackedParameter<double>( "maxX", 18.0 );
53  m_minY = iPset.getUntrackedParameter<double>( "minY", -18.0 );
54  m_maxY = iPset.getUntrackedParameter<double>( "maxY", 18.0 );
55  m_minZ = iPset.getUntrackedParameter<double>( "minZ", -18.0 );
56  m_maxZ = iPset.getUntrackedParameter<double>( "maxZ", 18.0 );
57 
58  m_xBaseDir = ( m_maxX - m_minX ) / m_mapDensityX;
59  m_yBaseDir = ( m_maxY - m_minY ) / m_mapDensityY;
60  m_zBaseDir = ( m_maxZ - m_minZ ) / m_mapDensityZ;
61 }
62 
64 {}
65 
66 void
68 {
69  iSetup.get<IdealMagneticFieldRecord>().get( m_mf );
70  m_mf.isValid() ? m_valid = true : m_valid = false;
71 
72  for( unsigned i = 0; i <= m_mapDensityX; ++i )
73  {
74  for( unsigned j = 0; j <= m_mapDensityY; ++j )
75  {
76  for( unsigned k = 0; k <= m_mapDensityZ; ++k )
77  {
78  // Prepare field position and get value.
79  double x = m_minX + m_xBaseDir * i;
80  double y = m_minY + m_yBaseDir * j;
81  double z = m_minZ + m_zBaseDir * k;
82  double pt[3] = { x, y, z };
83  double val[3];
84  evaluate( pt, val );
85  std::cout << "(" << x << ", " << y << ", " << z << ") " << val[0] << ":" << val[1] << ":" << val[2] << "; ";
86  }
87  std::cout << std::endl;
88  }
89  std::cout << std::endl;
90  }
91 }
92 
93 void
94 MFAnalyzer::evaluate (const double point [3], double field [3]) const
95 {
96  GlobalVector b = m_mf->inTesla( GlobalPoint( point[0], point[1], point[2] ));
97 
98  field [0] = b.x();
99  field [1] = b.y();
100  field [2] = b.z();
101 }
102 
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
double m_minX
Definition: MFAnalyzer.cc:32
void evaluate(const double point[3], double field[3]) const
Definition: MFAnalyzer.cc:94
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
MFAnalyzer(const edm::ParameterSet &)
Definition: MFAnalyzer.cc:45
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
T y() const
Definition: PV3DBase.h:57
double m_minZ
Definition: MFAnalyzer.cc:36
edm::ESHandle< MagneticField > m_mf
Definition: MFAnalyzer.cc:42
~MFAnalyzer(void)
Definition: MFAnalyzer.cc:63
double m_xBaseDir
Definition: MFAnalyzer.cc:38
int iEvent
Definition: GenABIO.cc:243
double m_maxX
Definition: MFAnalyzer.cc:33
double m_minY
Definition: MFAnalyzer.cc:34
Definition: DDAxes.h:10
unsigned m_mapDensityY
Definition: MFAnalyzer.cc:30
T z() const
Definition: PV3DBase.h:58
int j
Definition: DBlmapReader.cc:9
double m_maxZ
Definition: MFAnalyzer.cc:37
bool m_valid
Definition: MFAnalyzer.cc:41
int k[5][pyjets_maxn]
double m_zBaseDir
Definition: MFAnalyzer.cc:40
unsigned m_mapDensityX
Definition: MFAnalyzer.cc:29
double m_maxY
Definition: MFAnalyzer.cc:35
const T & get() const
Definition: EventSetup.h:55
double b
Definition: hdecay.h:120
double m_yBaseDir
Definition: MFAnalyzer.cc:39
tuple cout
Definition: gather_cfg.py:41
unsigned m_mapDensityZ
Definition: MFAnalyzer.cc:31
bool isValid() const
Definition: ESHandle.h:37
T x() const
Definition: PV3DBase.h:56
virtual void analyze(const edm::Event &, const edm::EventSetup &)
Definition: MFAnalyzer.cc:67
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5
bool valid(void)
Definition: MFAnalyzer.cc:24