CMS 3D CMS Logo
DisplayGeom.cc
Go to the documentation of this file.
1 
2 #include "FWCore/Framework/interface/Frameworkfwd.h"
3 #include "FWCore/Framework/interface/EDAnalyzer.h"
4 #include "FWCore/Framework/interface/MakerMacros.h"
5 
6 #include "FWCore/ParameterSet/interface/ParameterSet.h"
7 #include "FWCore/Utilities/interface/InputTag.h"
8 #include "FWCore/Framework/interface/Event.h"
9 #include "FWCore/Framework/interface/EventSetup.h"
10 #include "FWCore/Framework/interface/ESHandle.h"
11 #include "FWCore/Framework/interface/ESWatcher.h"
12 
13 #include "FWCore/ServiceRegistry/interface/Service.h"
14 
15 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
16 #include "MagneticField/Engine/interface/MagneticField.h"
17 
18 #include "Fireworks/Eve/interface/EveService.h"
19 
20 #include "Fireworks/Geometry/interface/DisplayGeomRecord.h"
21 
22 #include "DataFormats/TrackReco/interface/Track.h"
23 
24 #include "Fireworks/Tracks/interface/TrackUtils.h"
25 
26 #include <boost/algorithm/string/case_conv.hpp>
27 
28 #include "TEveManager.h"
29 #include "TEveTrack.h"
30 #include "TEveTrackPropagator.h"
31 
32 #include "TGeoManager.h"
33 #include "TGeoMatrix.h"
34 #include "TEveGeoNode.h"
35 #include "TEveTrans.h"
36 #include "TEveScene.h"
37 #include "TEveViewer.h"
38 #include "TVector.h"
39 #include "TGLScenePad.h"
40 #include "TGLRnrCtx.h"
41 #include "TEvePointSet.h"
42 #include "TRandom.h"
43 #include "TEveUtil.h"
44 
45 #include "TEveQuadSet.h"
46 #include "TEveStraightLineSet.h"
47 #include "TEveRGBAPalette.h"
48 #include "TSystem.h"
49 #include "TStyle.h"
50 // class decleration
51 //
52 
53 class DisplayGeom : public edm::EDAnalyzer {
54 public:
55  explicit DisplayGeom(const edm::ParameterSet&);
56  ~DisplayGeom() override;
57 
58 protected:
59  TEveGeoTopNode* make_node(const TString& path, Int_t vis_level, Bool_t global_cs);
60 
61 private:
62  void beginJob() override;
63  void analyze(const edm::Event&, const edm::EventSetup&) override;
64 
65  void endJob() override;
66 
67  edm::Service<EveService> m_eve;
68 
69  TEveElement* m_geomList;
70 
71  int m_level;
72 
73  bool m_MF;
74  int m_MF_component;
75  std::vector<double> m_MF_plane_d0;
76  std::vector<double> m_MF_plane_d1;
77  std::vector<double> m_MF_plane_d2;
78  int m_MF_plane_N1;
79  int m_MF_plane_N2;
80  int m_MF_plane_draw_dir;
81  bool m_MF_isPickable;
82 
83  edm::ESWatcher<DisplayGeomRecord> m_geomWatcher;
84 
85  void remakeGeometry(const DisplayGeomRecord& dgRec);
86 };
87 
88 DEFINE_FWK_MODULE(DisplayGeom);
89 
90 DisplayGeom::DisplayGeom(const edm::ParameterSet& iConfig)
91  : m_eve(), m_geomList(nullptr), m_MF_component(0), m_geomWatcher(this, &DisplayGeom::remakeGeometry) {
92  m_level = iConfig.getUntrackedParameter<int>("level", 2); // Geometry level to visualize
93 
94  m_MF = iConfig.getUntrackedParameter<int>("MF", false); // Show the MF geometry, instead of detector geometry
95 
96  std::string component = iConfig.getUntrackedParameter<std::string>("MF_component", "NONE");
97  boost::algorithm::to_upper(component);
98 
99  if (component == "NONE") {
100  m_MF_component = -1;
101  } else if (component == "ABSBZ") {
102  m_MF_component = 1;
103  } else if (component == "ABSBR") {
104  m_MF_component = 2;
105  } else if (component == "ABSBPHI") {
106  m_MF_component = 3;
107  } else if (component == "BR") {
108  m_MF_component = 4;
109  } else if (component == "BPHI") {
110  m_MF_component = 5;
111  } else { // Anything else -> |B|
112  m_MF_component = 0;
113  }
114 
115  if (m_MF_component == 0) {
116  m_MF_plane_d0 = iConfig.getUntrackedParameter<std::vector<double> >("MF_plane_d0", std::vector<double>(3, 0.0));
117  m_MF_plane_d1 = iConfig.getParameter<std::vector<double> >("MF_plane_d1");
118  m_MF_plane_d2 = iConfig.getParameter<std::vector<double> >("MF_plane_d2");
119 
120  m_MF_plane_N1 = iConfig.getUntrackedParameter<UInt_t>("MF_plane_N", 100);
121  m_MF_plane_N2 = iConfig.getUntrackedParameter<UInt_t>("MF_plane_N2", m_MF_plane_N1);
122 
123  m_MF_plane_draw_dir = iConfig.getUntrackedParameter<int>("MF_plane_draw_dir", true);
124  m_MF_isPickable = iConfig.getUntrackedParameter<bool>("MF_pickable", true);
125  }
126 }
127 
128 DisplayGeom::~DisplayGeom() {}
129 
130 //==============================================================================
131 // Protected helpers
132 //==============================================================================
133 
134 TEveGeoTopNode* DisplayGeom::make_node(const TString& path, Int_t vis_level, Bool_t global_cs) {
135  if (!gGeoManager->cd(path)) {
136  Warning("make_node", "Path '%s' not found.", path.Data());
137  return nullptr;
138  }
139 
140  TEveGeoTopNode* tn = new TEveGeoTopNode(gGeoManager, gGeoManager->GetCurrentNode());
141  tn->SetVisLevel(vis_level);
142  if (global_cs) {
143  tn->RefMainTrans().SetFrom(*gGeoManager->GetCurrentMatrix());
144  }
145  m_geomList->AddElement(tn);
146  gEve->AddToListTree(tn, true);
147  return tn;
148 }
149 
150 //==============================================================================
151 // member functions
152 //==============================================================================
153 
154 void DisplayGeom::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
155  if (m_eve) {
156  // Remake geometry if it has changed.
157  m_geomWatcher.check(iSetup);
158 
159  if (m_MF_component != -1) {
160  edm::ESHandle<MagneticField> field;
161  iSetup.get<IdealMagneticFieldRecord>().get(field);
162 
163  gStyle->SetPalette(1, nullptr);
164 
165  int minval = 0;
166  int maxval = 4000;
167  if (m_MF_component == 1) { //AbsBZ
168  minval = 0, maxval = 4000;
169  } else if (m_MF_component == 2) { //AbsBR
170  minval = 0, maxval = 4000;
171  } else if (m_MF_component == 3) { //AbsBphi
172  minval = 0, maxval = 1000;
173  } else if (m_MF_component == 4) { //BR
174  minval = -4000, maxval = 4000;
175  } else if (m_MF_component == 5) { //Bphi
176  minval = -1000, maxval = 1000;
177  }
178 
179  TEveRGBAPalette* pal = new TEveRGBAPalette(minval, maxval);
180 
181  TEveStraightLineSet* ls = nullptr;
182  if (m_MF_plane_draw_dir) {
183  new TEveStraightLineSet("MF_line_direction");
184  ls->SetPickable(false);
185  ls->SetLineColor(kGreen);
186  ls->SetMarkerColor(kGreen);
187  ls->SetMarkerStyle(1);
188  }
189 
190  TEveQuadSet* q = new TEveQuadSet("MF_quad_values");
191  q->Reset(TEveQuadSet::kQT_RectangleXY, kFALSE, 32);
192  q->SetOwnIds(kTRUE);
193  q->SetAlwaysSecSelect(true);
194  q->SetPickable(m_MF_isPickable);
195  q->SetPalette(pal);
196 
197  TEveVectorD v0(m_MF_plane_d0[0], m_MF_plane_d0[1], m_MF_plane_d0[2]);
198  TEveVectorD v01(m_MF_plane_d1[0], m_MF_plane_d1[1], m_MF_plane_d1[2]);
199  TEveVectorD v02(m_MF_plane_d2[0], m_MF_plane_d2[1], m_MF_plane_d2[2]);
200 
201  TEveVectorD b01 = (v01 - v0);
202  TEveVectorD b02 = (v02 - v0);
203  TEveVectorD b03 = b01.Cross(b02);
204 
205  TEveTrans trans;
206  trans.SetBaseVec(1, b01.fX, b01.fY, b01.fZ);
207  trans.SetBaseVec(2, b02.fX, b02.fY, b02.fZ);
208  trans.SetBaseVec(3, b03.fX, b03.fY, b03.fZ);
209  trans.SetPos(v0.Arr());
210  trans.OrtoNorm3();
211  q->SetTransMatrix(trans.Array());
212 
213  double w_step = b01.Mag() / m_MF_plane_N1;
214  double h_step = b02.Mag() / m_MF_plane_N2;
215 
216  q->SetDefWidth(w_step);
217  q->SetDefHeight(h_step);
218  TEveVectorD d1;
219  trans.GetBaseVec(1).GetXYZ(d1);
220  d1 *= w_step;
221  TEveVectorD d2;
222  trans.GetBaseVec(2).GetXYZ(d2);
223  d2 *= h_step;
224 
225  //d1.Print();
226  d2.Dump();
227  double line_step_size = TMath::Min(w_step, h_step);
228 
229  for (int i = 0; i < m_MF_plane_N1; i++) {
230  for (int j = 0; j < m_MF_plane_N2; j++) {
231  TEveVectorD p = d1 * Double_t(i) + d2 * Double_t(j) + v0;
232  GlobalPoint pos(p.fX, p.fY, p.fZ);
233  GlobalVector b = field->inTesla(pos) * 1000.; // in mT
234  float value = 0.;
235  if (m_MF_component == 0) { //BMOD
236  value = b.mag();
237  } else if (m_MF_component == 1) { //BZ
238  value = fabs(b.z());
239  } else if (m_MF_component == 2) { //ABSBR
240  value = fabs(b.x() * cos(pos.phi()) + b.y() * sin(pos.phi()));
241  } else if (m_MF_component == 3) { //ABSBPHI
242  value = fabs(-b.x() * sin(pos.phi()) + b.y() * cos(pos.phi()));
243  } else if (m_MF_component == 2) { //BR
244  value = b.x() * cos(pos.phi()) + b.y() * sin(pos.phi());
245  } else if (m_MF_component == 5) { //BPHI
246  value = -b.x() * sin(pos.phi()) + b.y() * cos(pos.phi());
247  }
248 
249  q->AddQuad(w_step * i, h_step * j);
250  q->QuadValue(value);
251  if (m_MF_isPickable)
252  q->QuadId(new TNamed(Form("Mag (%f, %f, %f) val = %f", b.x(), b.y(), b.z(), b.mag()), "Dong!"));
253 
254  if (ls) {
255  if (b.mag() > 1e-6) {
256  b.unit();
257  b *= line_step_size;
258  ls->AddLine(p.fX, p.fY, p.fZ, p.fX + b.x(), p.fY + b.y(), p.fZ + b.z());
259  } else {
260  ls->AddLine(p.fX, p.fY, p.fZ, p.fX + b.x(), p.fY + b.y(), p.fZ + b.z());
261  }
262 
263  ls->AddMarker(ls->GetLinePlex().Size() - 1, 0);
264  ls->AddMarker(i * m_MF_plane_N1 + j, 0);
265  }
266  }
267  }
268 
269  TEveScene* eps = gEve->SpawnNewScene("FillStyleScene");
270  gEve->GetDefaultViewer()->AddScene(eps);
271  eps->GetGLScene()->SetStyle(TGLRnrCtx::kFill);
272  eps->AddElement(q);
273  if (ls)
274  m_eve->AddElement(ls);
275  } else {
276  // // Add a test obj
277  // if (!gRandom)
278  // gRandom = new TRandom(0);
279  // TRandom& r= *gRandom;
280 
281  // Float_t s = 100;
282 
283  // TEvePointSet* ps = new TEvePointSet();
284  // ps->SetOwnIds(kTRUE);
285  // for(Int_t i = 0; i< 100; i++)
286  // {
287  // ps->SetNextPoint(r.Uniform(-s,s), r.Uniform(-s,s), r.Uniform(-s,s));
288  // ps->SetPointId(new TNamed(Form("Point %d", i), ""));
289  // }
290 
291  // ps->SetMarkerColor(TMath::Nint(r.Uniform(2, 9)));
292  // ps->SetMarkerSize(r.Uniform(1, 2));
293  // ps->SetMarkerStyle(4);
294  // m_eve->AddElement(ps);
295  }
296  }
297 }
298 
299 // ------------ method called once each job just before starting event loop ------------
300 void DisplayGeom::beginJob() {
301  if (m_eve) {
302  m_geomList = new TEveElementList("Display Geom");
303  m_eve->AddGlobalElement(m_geomList);
304  // m_eve->getManager()->GetGlobalScene()->GetGLScene()->SetStyle(TGLRnrCtx::kWireFrame);
305  }
306 }
307 
308 // ------------ method called once each job just after ending the event loop ------------
309 void DisplayGeom::endJob() {}
310 
311 //------------------------------------------------------------------------------
312 void DisplayGeom::remakeGeometry(const DisplayGeomRecord& dgRec) {
313  m_geomList->DestroyElements();
314 
315  edm::ESHandle<TGeoManager> geom;
316  dgRec.get(geom);
317  TEveGeoManagerHolder _tgeo(const_cast<TGeoManager*>(geom.product()));
318 
319  // To have a full one, all detectors in one top-node:
320  // make_node("/cms:World_1/cms:CMSE_1", 4, kTRUE);
321 
322  if (m_MF) {
323  make_node("/cms:World_1", m_level, kTRUE);
324  } else {
325  make_node("/cms:World_1/cms:CMSE_1/tracker:Tracker_1", m_level, kTRUE);
326  make_node("/cms:World_1/cms:CMSE_1/caloBase:CALO_1", m_level, kTRUE);
327  make_node("/cms:World_1/cms:CMSE_1/muonBase:MUON_1", m_level, kTRUE);
328  }
329 }
Vector3DBase
Definition: Vector3DBase.h:8
edm::ESWatcher::check
bool check(const edm::EventSetup &iSetup)
Definition: ESWatcher.h:52
eostools.ls
def ls(path, rec=False)
Definition: eostools.py:349
DisplayGeom::m_MF_plane_N2
int m_MF_plane_N2
Definition: DisplayGeom.cc:79
MagneticField::inTesla
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
V0Monitor_cff.v0
v0
Definition: V0Monitor_cff.py:7
mps_fire.i
i
Definition: mps_fire.py:355
DisplayGeom::m_MF_component
int m_MF_component
Definition: DisplayGeom.cc:74
DisplayGeom::m_MF_plane_draw_dir
int m_MF_plane_draw_dir
Definition: DisplayGeom.cc:80
edm::ESWatcher< DisplayGeomRecord >
ESHandle.h
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
DisplayGeom::remakeGeometry
void remakeGeometry(const DisplayGeomRecord &dgRec)
Definition: DisplayGeom.cc:312
pos
Definition: PixelAliasList.h:18
EveService::AddElement
void AddElement(TEveElement *el)
Definition: EveService.cc:261
data-class-funcs.q
q
Definition: data-class-funcs.py:169
DisplayGeom::m_MF_plane_d1
std::vector< double > m_MF_plane_d1
Definition: DisplayGeom.cc:76
DisplayGeom::m_eve
edm::Service< EveService > m_eve
Definition: DisplayGeom.cc:67
DisplayGeom::m_MF_isPickable
bool m_MF_isPickable
Definition: DisplayGeom.cc:81
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
EDAnalyzer.h
DisplayGeom::m_level
int m_level
Definition: DisplayGeom.cc:71
TrackUtils.h
DisplayGeom::DisplayGeom
DisplayGeom(const edm::ParameterSet &)
Definition: DisplayGeom.cc:90
edm::EDAnalyzer
Definition: EDAnalyzer.h:29
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
EveService::AddGlobalElement
void AddGlobalElement(TEveElement *el)
Definition: EveService.cc:263
IdealMagneticFieldRecord
Definition: IdealMagneticFieldRecord.h:11
MakerMacros.h
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Track.h
edm::EventSetup::get
T get() const
Definition: EventSetup.h:73
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
Service.h
DisplayGeom::make_node
TEveGeoTopNode * make_node(const TString &path, Int_t vis_level, Bool_t global_cs)
Definition: DisplayGeom.cc:134
IdealMagneticFieldRecord.h
edm::ESHandle< MagneticField >
relativeConstraints.geom
geom
Definition: relativeConstraints.py:72
Point3DBase< float, GlobalTag >
edm::eventsetup::DependentRecordImplementation::get
ProductT const & get(ESGetToken< ProductT, DepRecordT > const &iToken) const
Definition: DependentRecordImplementation.h:112
b
double b
Definition: hdecay.h:118
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
DisplayGeom
Definition: DisplayGeom.cc:53
DisplayGeom::m_MF_plane_d0
std::vector< double > m_MF_plane_d0
Definition: DisplayGeom.cc:75
DisplayGeom::endJob
void endJob() override
Definition: DisplayGeom.cc:309
edm::ParameterSet
Definition: ParameterSet.h:36
Event.h
DisplayGeom::m_MF_plane_d2
std::vector< double > m_MF_plane_d2
Definition: DisplayGeom.cc:77
edm::Service< EveService >
iEvent
int iEvent
Definition: GenABIO.cc:224
value
Definition: value.py:1
DisplayGeom::~DisplayGeom
~DisplayGeom() override
Definition: DisplayGeom.cc:128
MagneticField.h
edm::EventSetup
Definition: EventSetup.h:57
DisplayGeom::m_MF
bool m_MF
Definition: DisplayGeom.cc:73
get
#define get
EveService.h
InputTag.h
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
DisplayGeom::m_MF_plane_N1
int m_MF_plane_N1
Definition: DisplayGeom.cc:78
DisplayGeomRecord.h
Frameworkfwd.h
ESWatcher.h
DisplayGeom::m_geomList
TEveElement * m_geomList
Definition: DisplayGeom.cc:69
EventSetup.h
DisplayGeom::analyze
void analyze(const edm::Event &, const edm::EventSetup &) override
Definition: DisplayGeom.cc:154
castor_dqm_sourceclient_file_cfg.path
path
Definition: castor_dqm_sourceclient_file_cfg.py:37
ParameterSet.h
Min
T Min(T a, T b)
Definition: MathUtil.h:39
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
edm::Event
Definition: Event.h:73
d1
static constexpr float d1
Definition: L1EGammaCrystalsEmulatorProducer.cc:84
DisplayGeom::beginJob
void beginJob() override
Definition: DisplayGeom.cc:300
DisplayGeom::m_geomWatcher
edm::ESWatcher< DisplayGeomRecord > m_geomWatcher
Definition: DisplayGeom.cc:83
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
DisplayGeomRecord
Definition: DisplayGeomRecord.h:24
Generated for CMSSW Reference Manual by doxygen 1.8.16