Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
TMarker3DBox.cxx
Go to the documentation of this file.
1 // @(#)root/g3d:$Id$
2 // Author: Rene Brun , Olivier Couet 31/10/97
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2004, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #include "Riostream.h"
13 #include "TROOT.h"
14 #include "TView.h"
15 #include "TMarker3DBox.h"
16 #include "TVirtualPad.h"
17 #include "TH1.h"
18 #include "TH3.h"
19 #include "TBuffer3D.h"
20 #include "TBuffer3DTypes.h"
21 #include "TVirtualViewer3D.h"
22 #include "TGeometry.h"
23 #include "TClass.h"
24 #include "TMath.h"
25 
26 #include <assert.h>
27 
28 ClassImp(TMarker3DBox);
29 
30 /** \class TMarker3DBox
31 \ingroup g3d
32 A special 3-D marker designed for event display.
33 
34 It has the following parameters:
35  - fX: X coordinate of the center of the box
36  - fY: Y coordinate of the center of the box
37  - fZ: Z coordinate of the center of the box
38  - fDx: half length in X
39  - fDy: half length in Y
40  - fDz: half length in Z
41  - fTheta: Angle of box z axis with respect to main Z axis
42  - fPhi: Angle of box x axis with respect to main Xaxis
43  - fRefObject: A reference to an object
44 */
45 
46 ////////////////////////////////////////////////////////////////////////////////
47 /// Marker3DBox default constructor
48 
49 TMarker3DBox::TMarker3DBox()
50 {
51  fRefObject = 0;
52  fDx = 1;
53  fDy = 1;
54  fDz = 1;
55  fX = 0;
56  fY = 0;
57  fZ = 0;
58  fTheta = 0;
59  fPhi = 0;
60  SetBit(kTemporary,kFALSE);
61 }
62 
63 ////////////////////////////////////////////////////////////////////////////////
64 /// Marker3DBox normal constructor
65 
66 TMarker3DBox::TMarker3DBox( Float_t x, Float_t y, Float_t z,
67  Float_t dx, Float_t dy, Float_t dz,
68  Float_t theta, Float_t phi)
69  :TAttLine(1,1,1), TAttFill(1,0)
70 {
71  fDx = dx;
72  fDy = dy;
73  fDz = dz;
74  fX = x;
75  fY = y;
76  fZ = z;
77  fTheta = theta;
78  fPhi = phi;
79  fRefObject = 0;
80  SetBit(kTemporary,kFALSE);
81 }
82 
83 ////////////////////////////////////////////////////////////////////////////////
84 /// copy constructor
85 
86 TMarker3DBox::TMarker3DBox(const TMarker3DBox& m3d) :
87  TObject(m3d),
88  TAttLine(m3d),
89  TAttFill(m3d),
90  TAtt3D(m3d),
91  fX(m3d.fX),
92  fY(m3d.fY),
93  fZ(m3d.fZ),
94  fDx(m3d.fDx),
95  fDy(m3d.fDy),
96  fDz(m3d.fDz),
97  fTheta(m3d.fTheta),
98  fPhi(m3d.fPhi),
99  fRefObject(m3d.fRefObject)
100 {
101 }
102 
103 ////////////////////////////////////////////////////////////////////////////////
104 /// assignment operator
105 
106 TMarker3DBox& TMarker3DBox::operator=(const TMarker3DBox& m3d)
107 {
108  if(this!=&m3d) {
109  TObject::operator=(m3d);
110  TAttLine::operator=(m3d);
111  TAttFill::operator=(m3d);
112  TAtt3D::operator=(m3d);
113  fX=m3d.fX;
114  fY=m3d.fY;
115  fZ=m3d.fZ;
116  fDx=m3d.fDx;
117  fDy=m3d.fDy;
118  fDz=m3d.fDz;
119  fTheta=m3d.fTheta;
120  fPhi=m3d.fPhi;
121  fRefObject=m3d.fRefObject;
122  }
123  return *this;
124 }
125 
126 ////////////////////////////////////////////////////////////////////////////////
127 /// Marker3DBox shape default destructor
128 
129 TMarker3DBox::~TMarker3DBox()
130 {
131 }
132 
133 
134 ////////////////////////////////////////////////////////////////////////////////
135 /// Compute distance from point px,py to a Marker3DBox
136 ///
137 /// Compute the closest distance of approach from point px,py to each corner
138 /// point of the Marker3DBox.
139 
140 Int_t TMarker3DBox::DistancetoPrimitive(Int_t px, Int_t py)
141 {
142  const Int_t numPoints = 8;
143  Int_t dist = 9999;
144  Double_t points[3*numPoints];
145 
146  TView *view = gPad->GetView();
147  if (!view) return dist;
148  const Int_t seg1[12] = {0,1,2,3,4,5,6,7,0,1,2,3};
149  const Int_t seg2[12] = {1,2,3,0,5,6,7,4,4,5,6,7};
150 
151  SetPoints(points);
152 
153  Int_t i, i1, i2, dsegment;
154  Double_t x1,y1,x2,y2;
155  Double_t xndc[3];
156  for (i = 0; i < 12; i++) {
157  i1 = 3*seg1[i];
158  view->WCtoNDC(&points[i1], xndc);
159  x1 = xndc[0];
160  y1 = xndc[1];
161 
162  i2 = 3*seg2[i];
163  view->WCtoNDC(&points[i2], xndc);
164  x2 = xndc[0];
165  y2 = xndc[1];
166  dsegment = DistancetoLine(px,py,x1,y1,x2,y2);
167  if (dsegment < dist) dist = dsegment;
168  }
169  if (dist < 5) {
170  gPad->SetCursor(kCross);
171  if (fRefObject) {gPad->SetSelected(fRefObject); return 0;}
172  }
173  return dist;
174 }
175 
176 ////////////////////////////////////////////////////////////////////////////////
177 /// Execute action corresponding to one event
178 ///
179 /// This member function must be implemented to realize the action
180 /// corresponding to the mouse click on the object in the window
181 
182 void TMarker3DBox::ExecuteEvent(Int_t event, Int_t px, Int_t py)
183 {
184  if (!gPad) return;
185  if (gPad->GetView()) gPad->GetView()->ExecuteRotateView(event, px, py);
186 }
187 
188 ////////////////////////////////////////////////////////////////////////////////
189 /// Paint marker 3D box.
190 
191 void TMarker3DBox::Paint(Option_t * /* option */ )
192 {
193  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
194 
195  buffer.ClearSectionsValid();
196 
197  // Section kCore
198 
199  // If we are just a temporary object then no 'real object' to
200  // pass to viewer
201  if (TestBit(kTemporary)) {
202  buffer.fID = 0;
203  } else {
204  buffer.fID = this;
205  }
206  buffer.fColor = GetLineColor();
207  buffer.fTransparency = 0;
208  buffer.fLocalFrame = kFALSE;
209  buffer.SetSectionsValid(TBuffer3D::kCore);
210 
211  // We fill kCore and kRawSizes on first pass and try with viewer
212  TVirtualViewer3D * viewer3D = gPad->GetViewer3D();
213  if (!viewer3D) return;
214  Int_t reqSections = viewer3D->AddObject(buffer);
215  if (reqSections == TBuffer3D::kNone) {
216  return;
217  }
218 
219  if (reqSections & TBuffer3D::kRawSizes) {
220  Int_t nbPnts = 8;
221  Int_t nbSegs = 12;
222  Int_t nbPols = 6;
223  if (!buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
224  return;
225  }
226  buffer.SetSectionsValid(TBuffer3D::kRawSizes);
227  }
228 
229  if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
230  // Points
231  SetPoints(buffer.fPnts);
232 
233  // Transform points
234  if (gGeometry && !buffer.fLocalFrame) {
235  Double_t dlocal[3];
236  Double_t dmaster[3];
237  for (UInt_t j=0; j<buffer.NbPnts(); j++) {
238  dlocal[0] = buffer.fPnts[3*j];
239  dlocal[1] = buffer.fPnts[3*j+1];
240  dlocal[2] = buffer.fPnts[3*j+2];
241  gGeometry->Local2Master(&dlocal[0],&dmaster[0]);
242  buffer.fPnts[3*j] = dmaster[0];
243  buffer.fPnts[3*j+1] = dmaster[1];
244  buffer.fPnts[3*j+2] = dmaster[2];
245  }
246  }
247 
248  // Basic colors: 0, 1, ... 8
249  Int_t c = (((GetLineColor()) %8) -1) * 4;
250  if (c < 0) c = 0;
251 
252  // Segments
253  buffer.fSegs[ 0] = c ; buffer.fSegs[ 1] = 0 ; buffer.fSegs[ 2] = 1;
254  buffer.fSegs[ 3] = c+1 ; buffer.fSegs[ 4] = 1 ; buffer.fSegs[ 5] = 2;
255  buffer.fSegs[ 6] = c+1 ; buffer.fSegs[ 7] = 2 ; buffer.fSegs[ 8] = 3;
256  buffer.fSegs[ 9] = c ; buffer.fSegs[10] = 3 ; buffer.fSegs[11] = 0;
257  buffer.fSegs[12] = c+2 ; buffer.fSegs[13] = 4 ; buffer.fSegs[14] = 5;
258  buffer.fSegs[15] = c+2 ; buffer.fSegs[16] = 5 ; buffer.fSegs[17] = 6;
259  buffer.fSegs[18] = c+3 ; buffer.fSegs[19] = 6 ; buffer.fSegs[20] = 7;
260  buffer.fSegs[21] = c+3 ; buffer.fSegs[22] = 7 ; buffer.fSegs[23] = 4;
261  buffer.fSegs[24] = c ; buffer.fSegs[25] = 0 ; buffer.fSegs[26] = 4;
262  buffer.fSegs[27] = c+2 ; buffer.fSegs[28] = 1 ; buffer.fSegs[29] = 5;
263  buffer.fSegs[30] = c+1 ; buffer.fSegs[31] = 2 ; buffer.fSegs[32] = 6;
264  buffer.fSegs[33] = c+3 ; buffer.fSegs[34] = 3 ; buffer.fSegs[35] = 7;
265 
266  // Polygons
267  buffer.fPols[ 0] = c ; buffer.fPols[ 1] = 4 ; buffer.fPols[ 2] = 0;
268  buffer.fPols[ 3] = 9 ; buffer.fPols[ 4] = 4 ; buffer.fPols[ 5] = 8;
269  buffer.fPols[ 6] = c+1 ; buffer.fPols[ 7] = 4 ; buffer.fPols[ 8] = 1;
270  buffer.fPols[ 9] = 10 ; buffer.fPols[10] = 5 ; buffer.fPols[11] = 9;
271  buffer.fPols[12] = c ; buffer.fPols[13] = 4 ; buffer.fPols[14] = 2;
272  buffer.fPols[15] = 11 ; buffer.fPols[16] = 6 ; buffer.fPols[17] = 10;
273  buffer.fPols[18] = c+1 ; buffer.fPols[19] = 4 ; buffer.fPols[20] = 3;
274  buffer.fPols[21] = 8 ; buffer.fPols[22] = 7 ; buffer.fPols[23] = 11;
275  buffer.fPols[24] = c+2 ; buffer.fPols[25] = 4 ; buffer.fPols[26] = 0;
276  buffer.fPols[27] = 3 ; buffer.fPols[28] = 2 ; buffer.fPols[29] = 1;
277  buffer.fPols[30] = c+3 ; buffer.fPols[31] = 4 ; buffer.fPols[32] = 4;
278  buffer.fPols[33] = 5 ; buffer.fPols[34] = 6 ; buffer.fPols[35] = 7;
279 
280  buffer.SetSectionsValid(TBuffer3D::kRaw);
281 
282  TAttLine::Modify();
283  TAttFill::Modify();
284  }
285 
286  viewer3D->AddObject(buffer);
287 }
288 
289 ////////////////////////////////////////////////////////////////////////////////
290 /// Paint 3-d histogram h with marker3dboxes
291 
292 void TMarker3DBox::PaintH3(TH1 *h, Option_t *option)
293 {
294  Int_t bin,ix,iy,iz;
295  Double_t xmin,xmax,ymin,ymax,zmin,zmax,wmin,wmax,w;
296  TAxis *xaxis = h->GetXaxis();
297  TAxis *yaxis = h->GetYaxis();
298  TAxis *zaxis = h->GetZaxis();
299 
300  wmin = h->GetMinimum();
301  wmax = h->GetMaximum();
302 
303  //Create or modify 3-d view object
304  TView *view = gPad->GetView();
305  if (!view) {
306  gPad->Range(-1,-1,1,1);
307  view = TView::CreateView(1,0,0);
308  if (!view) return;
309  }
310  view->SetRange(xaxis->GetBinLowEdge(xaxis->GetFirst()),
311  yaxis->GetBinLowEdge(yaxis->GetFirst()),
312  zaxis->GetBinLowEdge(zaxis->GetFirst()),
313  xaxis->GetBinUpEdge(xaxis->GetLast()),
314  yaxis->GetBinUpEdge(yaxis->GetLast()),
315  zaxis->GetBinUpEdge(zaxis->GetLast()));
316 
317  view->PadRange(gPad->GetFrameFillColor());
318 
319  //Draw TMarker3DBox with size proportional to cell content
320  TMarker3DBox m3;
321  m3.SetBit(kTemporary,kTRUE);
322  m3.SetRefObject(h);
323  m3.SetDirection(0,0);
324  m3.SetLineColor(h->GetMarkerColor());
325  Double_t scale;
326  for (ix=xaxis->GetFirst();ix<=xaxis->GetLast();ix++) {
327  xmin = h->GetXaxis()->GetBinLowEdge(ix);
328  xmax = xmin + h->GetXaxis()->GetBinWidth(ix);
329  for (iy=yaxis->GetFirst();iy<=yaxis->GetLast();iy++) {
330  ymin = h->GetYaxis()->GetBinLowEdge(iy);
331  ymax = ymin + h->GetYaxis()->GetBinWidth(iy);
332  for (iz=zaxis->GetFirst();iz<=zaxis->GetLast();iz++) {
333  zmin = h->GetZaxis()->GetBinLowEdge(iz);
334  zmax = zmin + h->GetZaxis()->GetBinWidth(iz);
335  bin = h->GetBin(ix,iy,iz);
336  w = h->GetBinContent(bin);
337  if (w < wmin) continue;
338  if (w > wmax) w = wmax;
339  scale = (TMath::Power((w-wmin)/(wmax-wmin),1./3.))/2.;
340  if (scale == 0) continue;
341  m3.SetPosition(0.5*(xmin+xmax),0.5*(ymin+ymax),0.5*(zmin+zmax));
342  m3.SetSize(scale*(xmax-xmin),scale*(ymax-ymin),scale*(zmax-zmin));
343  m3.Paint(option);
344  }
345  }
346  }
347 }
348 
349 ////////////////////////////////////////////////////////////////////////////////
350 /// Save primitive as a C++ statement(s) on output stream out
351 
352 void TMarker3DBox::SavePrimitive(std::ostream &out, Option_t * /*= ""*/)
353 {
354  out<<" "<<std::endl;
355  if (gROOT->ClassSaved(TMarker3DBox::Class())) {
356  out<<" ";
357  } else {
358  out<<" TMarker3DBox *";
359  }
360  out<<"marker3DBox = new TMarker3DBox("<<fX<<","
361  <<fY<<","
362  <<fZ<<","
363  <<fDx<<","
364  <<fDy<<","
365  <<fDz<<","
366  <<fTheta<<","
367  <<fPhi<<");"<<std::endl;
368 
369  SaveLineAttributes(out,"marker3DBox",1,1,1);
370  SaveFillAttributes(out,"marker3DBox",1,0);
371 
372  out<<" marker3DBox->Draw();"<<std::endl;
373 }
374 
375 ////////////////////////////////////////////////////////////////////////////////
376 /// Set direction.
377 
378 void TMarker3DBox::SetDirection(Float_t theta, Float_t phi)
379 {
380  fTheta = theta;
381  fPhi = phi;
382 }
383 
384 ////////////////////////////////////////////////////////////////////////////////
385 /// Set size.
386 
387 void TMarker3DBox::SetSize(Float_t dx, Float_t dy, Float_t dz)
388 {
389  fDx = dx;
390  fDy = dy;
391  fDz = dz;
392 }
393 
394 ////////////////////////////////////////////////////////////////////////////////
395 /// Set position.
396 
397 void TMarker3DBox::SetPosition(Float_t x, Float_t y, Float_t z)
398 {
399  fX = x;
400  fY = y;
401  fZ = z;
402 }
403 
404 ////////////////////////////////////////////////////////////////////////////////
405 /// Set points.
406 
407 void TMarker3DBox::SetPoints(Double_t *points) const
408 {
409  if (points) {
410  points[ 0] = -fDx ; points[ 1] = -fDy ; points[ 2] = -fDz;
411  points[ 3] = -fDx ; points[ 4] = fDy ; points[ 5] = -fDz;
412  points[ 6] = fDx ; points[ 7] = fDy ; points[ 8] = -fDz;
413  points[ 9] = fDx ; points[10] = -fDy ; points[11] = -fDz;
414  points[12] = -fDx ; points[13] = -fDy ; points[14] = fDz;
415  points[15] = -fDx ; points[16] = fDy ; points[17] = fDz;
416  points[18] = fDx ; points[19] = fDy ; points[20] = fDz;
417  points[21] = fDx ; points[22] = -fDy ; points[23] = fDz;
418 
419  Double_t x, y, z;
420  const Double_t kPI = TMath::Pi();
421  Double_t theta = fTheta*kPI/180;
422  Double_t phi = fPhi*kPI/180;
423  Double_t sinth = TMath::Sin(theta);
424  Double_t costh = TMath::Cos(theta);
425  Double_t sinfi = TMath::Sin(phi);
426  Double_t cosfi = TMath::Cos(phi);
427 
428  // Matrix to convert from fruit frame to master frame
429  Double_t m[9];
430  m[0] = costh * cosfi; m[1] = -sinfi; m[2] = sinth*cosfi;
431  m[3] = costh * sinfi; m[4] = cosfi; m[5] = sinth*sinfi;
432  m[6] = -sinth; m[7] = 0; m[8] = costh;
433  for (Int_t i = 0; i < 8; i++) {
434  x = points[3*i];
435  y = points[3*i+1];
436  z = points[3*i+2];
437 
438  points[3*i] = fX + m[0] * x + m[1] * y + m[2] * z;
439  points[3*i+1] = fY + m[3] * x + m[4] * y + m[5] * z;
440  points[3*i+2] = fZ + m[6] * x + m[7] * y + m[8] * z;
441  }
442  }
443 }
444 
445 ////////////////////////////////////////////////////////////////////////////////
446 /// Stream an object of class TMarker3DBox.
447 
448 void TMarker3DBox::Streamer(TBuffer &R__b)
449 {
450  if (R__b.IsReading()) {
451  UInt_t R__s, R__c;
452  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
453  if (R__v > 1) {
454  R__b.ReadClassBuffer(TMarker3DBox::Class(), this, R__v, R__s, R__c);
455  return;
456  }
457  //====process old versions before automatic schema evolution
458  TObject::Streamer(R__b);
459  TAttLine::Streamer(R__b);
460  TAttFill::Streamer(R__b);
461  TAtt3D::Streamer(R__b);
462  R__b >> fX;
463  R__b >> fY;
464  R__b >> fZ;
465  R__b >> fDx;
466  R__b >> fDy;
467  R__b >> fDz;
468  R__b >> fTheta;
469  R__b >> fPhi;
470  R__b >> fRefObject;
471  R__b.CheckByteCount(R__s, R__c, TMarker3DBox::IsA());
472  //====end of old versions
473 
474  } else {
475  R__b.WriteClassBuffer(TMarker3DBox::Class(),this);
476  }
477 }