CMS 3D CMS Logo

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