Logo ROOT  
Reference Guide
TView3D.cxx
Go to the documentation of this file.
1// @(#)root/g3d:$Id$
2// Author: Rene Brun, Nenad Buncic, Evgueni Tcherniaev, Olivier Couet 18/08/95
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, 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 "RConfigure.h"
13
14#include "TVirtualPad.h"
15#include "TView3D.h"
16#include "TAxis3D.h"
17#include "TPolyLine3D.h"
18#include "TVirtualX.h"
19#include "TROOT.h"
20#include "TBuffer.h"
21#include "TClass.h"
22#include "TList.h"
23#include "TPluginManager.h"
24#include "TMath.h"
25
26// Remove when TView3Der3DPad fix in ExecuteRotateView() is removed
27#include "TVirtualViewer3D.h"
28
30
31//const Int_t kPerspective = BIT(14);
32
33const Int_t kCARTESIAN = 1;
34const Int_t kPOLAR = 2;
35const Double_t kRad = 3.14159265358979323846/180.0;
36
37/** \class TView3D
38\ingroup g3d
39The 3D view class.
40
41This package was originally written by Evgueni Tcherniaev from IHEP/Protvino.
42
43The original Fortran implementation was adapted to HIGZ/PAW by Olivier Couet and
44Evgueni Tcherniaev.
45
46This View class is a subset of the original system. It has been converted to a
47C++ class by Rene Brun.
48
49TView3D creates a 3-D view in the current pad. In this 3D view Lego and Surface
50plots can be drawn and also 3D polyline and markers. Most of the time a TView3D
51is created automatically when a 3D object needs to be painted in a pad (for
52instance a Lego or a Surface plot).
53
54In some case a TView3D should be explicitly. For instance to paint a 3D simple
55scene composed of simple objects like polylines and polymarkers.
56The following macro gives an example:
57
58Begin_Macro(source)
59{
60 cV3D = new TCanvas("cV3D","PolyLine3D & PolyMarker3D Window",200,10,500,500);
61
62 // Creating a view
63 TView3D *view = (TView3D*) TView::CreateView(1);
64 view->SetRange(5,5,5,25,25,25);
65
66 // Create a first PolyLine3D
67 TPolyLine3D *pl3d1 = new TPolyLine3D(6);
68 pl3d1->SetPoint(0, 10, 20, 10);
69 pl3d1->SetPoint(1, 15, 15, 15);
70 pl3d1->SetPoint(2, 20, 20, 20);
71 pl3d1->SetPoint(3, 20, 10, 20);
72 pl3d1->SetPoint(4, 10, 10, 20);
73 pl3d1->SetPoint(5, 10, 10, 10);
74
75 // Create a first PolyMarker3D
76 TPolyMarker3D *pm3d1 = new TPolyMarker3D(9);
77 pm3d1->SetPoint( 0, 10, 10, 10);
78 pm3d1->SetPoint( 1, 20, 20, 20);
79 pm3d1->SetPoint( 2, 10, 20, 20);
80 pm3d1->SetPoint( 3, 10, 10, 20);
81 pm3d1->SetPoint( 4, 20, 20, 10);
82 pm3d1->SetPoint( 5, 20, 10, 10);
83 pm3d1->SetPoint( 6, 20, 10, 20);
84 pm3d1->SetPoint( 7, 10, 20, 10);
85 pm3d1->SetPoint( 8, 15, 15, 15);
86 pm3d1->SetMarkerSize(2);
87 pm3d1->SetMarkerColor(4);
88 pm3d1->SetMarkerStyle(2);
89
90 // Draw
91 pl3d1->Draw();
92 pm3d1->Draw();
93}
94End_Macro
95
96
97Several coordinate systems are available:
98
99 - Cartesian
100 - Polar
101 - Cylindrical
102 - Spherical
103 - PseudoRapidity/Phi
104*/
105
106////////////////////////////////////////////////////////////////////////////////
107/// Default constructor
108
110{
111 fSystem = 0;
112 fOutline = 0;
116
117 fPsi = 0;
118 Int_t i;
119 for (i = 0; i < 3; i++) {
120 fRmin[i] = 0;
121 fRmax[i] = 1;
122 fX1[i] = fX2[i] = fY1[i] = fY2[i] = fZ1[i] = fZ2[i] = 0;
123 }
124
125 if (gPad) {
126 fLongitude = -90 - gPad->GetPhi();
127 fLatitude = 90 - gPad->GetTheta();
128 } else {
129 fLongitude = 0;
130 fLatitude = 0;
131 }
132 Int_t irep = 1;
134}
135
136////////////////////////////////////////////////////////////////////////////////
137/// TView3D constructor
138///
139/// Creates a 3-D view in the current pad
140/// rmin[3], rmax[3] are the limits of the object depending on
141/// the selected coordinate system
142///
143/// Before drawing a 3-D object in a pad, a 3-D view must be created.
144/// Note that a view is automatically created when drawing legos or surfaces.
145///
146/// The coordinate system is selected via system:
147/// - system = 1 Cartesian
148/// - system = 2 Polar
149/// - system = 3 Cylindrical
150/// - system = 4 Spherical
151/// - system = 5 PseudoRapidity/Phi
152
153TView3D::TView3D(Int_t system, const Double_t *rmin, const Double_t *rmax) : TView()
154{
155 Int_t irep;
156
158
159 fSystem = system;
160 fOutline = 0;
164
165 if (system == kCARTESIAN || system == kPOLAR || system == 11) fPsi = 0;
166 else fPsi = 90;
167
168 // By default pad range in 3-D view is (-1,-1,1,1), so ...
169 if (gPad) gPad->Range(-1, -1, 1, 1);
171
172 Int_t i;
173 for (i = 0; i < 3; i++) {
174 if (rmin) fRmin[i] = rmin[i];
175 else fRmin[i] = 0;
176 if (rmax) fRmax[i] = rmax[i];
177 else fRmax[i] = 1;
178 fX1[i] = fX2[i] = fY1[i] = fY2[i] = fZ1[i] = fZ2[i] = 0;
179 }
180
181 if (gPad) {
182 fLongitude = -90 - gPad->GetPhi();
183 fLatitude = 90 - gPad->GetTheta();
184 } else {
185 fLongitude = 0;
186 fLatitude = 0;
187 }
189
190 if (gPad) gPad->SetView(this);
191 if (system == 11) SetPerspective();
192}
193
194////////////////////////////////////////////////////////////////////////////////
195/// Copy constructor.
196
198 :TView(tv),
199 fLatitude(tv.fLatitude),
200 fLongitude(tv.fLongitude),
201 fPsi(tv.fPsi),
202 fDview(tv.fDview),
203 fDproj(tv.fDproj),
204 fUpix(tv.fUpix),
205 fVpix(tv.fVpix),
206 fSystem(tv.fSystem),
207 fOutline(tv.fOutline),
208 fDefaultOutline(tv.fDefaultOutline),
209 fAutoRange(tv.fAutoRange),
210 fChanged(tv.fChanged)
211{
212 for (Int_t i=0; i<16; i++) {
213 fTN[i]=tv.fTN[i];
214 fTB[i]=tv.fTB[i];
215 fTnorm[i]=tv.fTnorm[i];
216 fTback[i]=tv.fTback[i];
217 }
218 for(Int_t i=0; i<3; i++) {
219 fRmax[i]=tv.fRmax[i];
220 fRmin[i]=tv.fRmin[i];
221 fX1[i]=tv.fX1[i];
222 fX2[i]=tv.fX2[i];
223 fY1[i]=tv.fY1[i];
224 fY2[i]=tv.fY2[i];
225 fZ1[i]=tv.fZ1[i];
226 fZ2[i]=tv.fZ2[i];
227 }
228 for(Int_t i=0; i<4; i++)
229 fUVcoord[i]=tv.fUVcoord[i];
230}
231
232////////////////////////////////////////////////////////////////////////////////
233/// Assignment operator.
234
236{
237 if (this!=&tv) {
241 fPsi=tv.fPsi;
242 fDview=tv.fDview;
243 fDproj=tv.fDproj;
244 fUpix=tv.fUpix;
245 fVpix=tv.fVpix;
246 fSystem=tv.fSystem;
251 for(Int_t i=0; i<16; i++) {
252 fTN[i]=tv.fTN[i];
253 fTB[i]=tv.fTB[i];
254 fTnorm[i]=tv.fTnorm[i];
255 fTback[i]=tv.fTback[i];
256 }
257 for(Int_t i=0; i<3; i++) {
258 fRmax[i]=tv.fRmax[i];
259 fRmin[i]=tv.fRmin[i];
260 fX1[i]=tv.fX1[i];
261 fX2[i]=tv.fX2[i];
262 fY1[i]=tv.fY1[i];
263 fY2[i]=tv.fY2[i];
264 fZ1[i]=tv.fZ1[i];
265 fZ2[i]=tv.fZ2[i];
266 }
267 for(Int_t i=0; i<4; i++)
268 fUVcoord[i]=tv.fUVcoord[i];
269 }
270 return *this;
271}
272
273////////////////////////////////////////////////////////////////////////////////
274/// TView3D default destructor.
275
277{
278 if (fOutline) fOutline->Delete();
279 delete fOutline;
280 fOutline = 0;
281}
282
283////////////////////////////////////////////////////////////////////////////////
284/// Define axis vertices.
285///
286/// Input:
287/// - ANG - angle between X and Y axis (not used anymore)
288///
289/// Output:
290/// - AV(3,8) - axis vertices
291/// - IX1 - 1st point of X-axis (x-min)
292/// - IX2 - 2nd point of X-axis (x-max)
293/// - IY1 - 1st point of Y-axis (y-min)
294/// - IY2 - 2nd point of Y-axis (y-max)
295/// - IZ1 - 1st point of Z-axis (z-min)
296/// - IZ2 - 2nd point of Z-axis (z-max)
297
298void TView3D::AxisVertex(Double_t, Double_t *av, Int_t &ix1, Int_t &ix2, Int_t &iy1, Int_t &iy2, Int_t &iz1, Int_t &iz2)
299{
300 Double_t p[8][3] = {
301 { fRmin[0], fRmin[1], fRmin[2] },
302 { fRmax[0], fRmin[1], fRmin[2] },
303 { fRmax[0], fRmax[1], fRmin[2] },
304 { fRmin[0], fRmax[1], fRmin[2] },
305 { fRmin[0], fRmin[1], fRmax[2] },
306 { fRmax[0], fRmin[1], fRmax[2] },
307 { fRmax[0], fRmax[1], fRmax[2] },
308 { fRmin[0], fRmax[1], fRmax[2] }
309 };
310 Int_t inodes[4][8] = {
311 { 2,3,4,1, 6,7,8,5 }, // x+, y+
312 { 3,4,1,2, 7,8,5,6 }, // x-, y+
313 { 1,2,3,4, 5,6,7,8 }, // x+, y-
314 { 4,1,2,3, 8,5,6,7 } // x-, y-
315 };
316 Int_t ixyminmax[16][4] = { // 8
317 { 3,2, 1,2 }, // x+, y+, z+, z-up 5 / \ 7
318 { 2,1, 3,2 }, // x-, y+, z+, z-up |\6/|
319 { 1,2, 2,3 }, // x+, y-, z+, z-up | | | Top view
320 { 2,3, 2,1 }, // x-, y-, z+, z-up 1 \|/ 3
321 // 2 6
322 { 4,1, 4,3 }, // x+, y+, z-, z-up 5 /|\ 7
323 { 3,4, 4,1 }, // x-, y+, z-, z-up | | |
324 { 4,3, 1,4 }, // x+, y-, z-, z-up Bottom view |/2\|
325 { 1,4, 3,4 }, // x-, y-, z-, z-up 1 \ / 3
326 // 2 4
327 { 8,5, 8,7 }, // x+, y+, z+, z-down 3 /|\ 1
328 { 7,8, 8,5 }, // x-, y+, z+, z-down | | |
329 { 8,7, 5,8 }, // x+, y-, z+, z-down |/6\| Bottom view
330 { 5,8, 7,8 }, // x-, y-, z+, z-down 7 \ / 5
331 // 8 4
332 { 7,6, 5,6 }, // x+, y+, z-, z-down 3 / \ 1
333 { 6,5, 7,6 }, // x-, y+, z-, z-down |\2/|
334 { 5,6, 6,7 }, // x+, y-, z-, z-down Top view | | |
335 { 6,7, 6,5 } // x-, y-, z-, z-down 7 \|/ 5
336 }; // 6
337
338 // Set vertices
339 Int_t icase = 0;
340 if (fTnorm[ 8] <= 0) icase += 1; // z projection of (1,0,0)
341 if (fTnorm[ 9] <= 0) icase += 2; // z projection of (0,1,0)
342 for (Int_t i=0; i<8; ++i) {
343 Int_t k = inodes[icase][i] - 1;
344 av[i*3+0] = p[k][0];
345 av[i*3+1] = p[k][1];
346 av[i*3+2] = p[k][2];
347 }
348
349 // Set indices for min and max
350 if (fTnorm[10] < 0) icase += 4; // z projection of (0,0,1)
351 if (fTnorm[ 6] < 0) icase += 8; // y projection of (0,0,1)
352 ix1 = ixyminmax[icase][0];
353 ix2 = ixyminmax[icase][1];
354 iy1 = ixyminmax[icase][2];
355 iy2 = ixyminmax[icase][3];
356 iz1 = (icase < 8) ? 1 : 3;
357 iz2 = (icase < 8) ? 5 : 7;
358}
359
360////////////////////////////////////////////////////////////////////////////////
361/// Define perspective view.
362///
363/// Compute transformation matrix from world coordinates
364/// to normalized coordinates (-1 to +1)
365///
366/// Input :
367/// - theta, phi - spherical angles giving the direction of projection
368/// - psi - screen rotation angle
369/// - cov[3] - center of view
370/// - dview - distance from COV to COP (center of projection)
371/// - umin, umax, vmin, vmax - view window in projection plane
372/// - dproj - distance from COP to projection plane
373/// - bcut, fcut - backward/forward range w.r.t projection plane (fcut<=0)
374///
375/// Output :
376/// - nper[16] - normalizing transformation
377///
378/// compute tr+rot to get COV in origin, view vector parallel to -Z axis, up
379/// vector parallel to Y.
380///
381/// ~~~ {.cpp}
382/// ^Yv UP ^ proj. plane
383/// | | /|
384/// | | / |
385/// | dproj / x--- center of window (COW)
386/// COV |----------|--x--|------------> Zv
387/// / | VRP'z
388/// / ---> | /
389/// / VPN |/
390/// Xv
391/// ~~~
392///
393/// 1. translate COP to origin of MARS : Tper = T(-copx, -copy, -copz)
394/// 2. rotate VPN : R = Rz(-psi)*Rx(-theta)*Rz(-phi) (inverse Euler)
395/// 3. left-handed screen reference to right-handed one of MARS : Trl
396///
397/// T12 = Tper*R*Trl
398
399
401{
402 Double_t t12[16];
403 Double_t cov[3];
404 Int_t i;
405 for (i=0; i<3; i++) cov[i] = 0.5*(fRmax[i]+fRmin[i]);
406
413
414 t12[0] = c1*c3 - s1*c2*s3;
415 t12[4] = c1*s3 + s1*c2*c3;
416 t12[8] = s1*s2;
417 t12[3] = 0;
418
419 t12[1] = -s1*c3 - c1*c2*s3;
420 t12[5] = -s1*s3 + c1*c2*c3;
421 t12[9] = c1*s2;
422 t12[7] = 0;
423
424 t12[2] = s2*s3;
425 t12[6] = -s2*c3;
426 t12[10] = c2; // contains Trl
427 t12[11] = 0;
428
429 // translate with -COP (before rotation):
430 t12[12] = -(cov[0]*t12[0]+cov[1]*t12[4]+cov[2]*t12[8]);
431 t12[13] = -(cov[0]*t12[1]+cov[1]*t12[5]+cov[2]*t12[9]);
432 t12[14] = -(cov[0]*t12[2]+cov[1]*t12[6]+cov[2]*t12[10]);
433 t12[15] = 1;
434
435 // translate with (0, 0, -dview) after rotation
436
437 t12[14] -= fDview;
438
439 // reflection on Z :
440 t12[2] *= -1;
441 t12[6] *= -1;
442 t12[10] *= -1;
443 t12[14] *= -1;
444
445 // Now we shear the center of window from (0.5*(umin+umax), 0.5*(vmin+vmax), dproj)
446 // to (0, 0, dproj)
447
448 Double_t a2 = -fUVcoord[0]/fDproj; // shear coef. on x
449 Double_t b2 = -fUVcoord[1]/fDproj; // shear coef. on y
450
451 // | 1 0 0 0 |
452 // SHz(a2,b2) = | 0 1 0 0 |
453 // | a2 b2 1 0 |
454 // | 0 0 0 1 |
455
456 fTnorm[0] = t12[0] + a2*t12[2];
457 fTnorm[1] = t12[1] + b2*t12[2];
458 fTnorm[2] = t12[2];
459 fTnorm[3] = 0;
460
461 fTnorm[4] = t12[4] + a2*t12[6];
462 fTnorm[5] = t12[5] + b2*t12[6];
463 fTnorm[6] = t12[6];
464 fTnorm[7] = 0;
465
466 fTnorm[8] = t12[8] + a2*t12[10];
467 fTnorm[9] = t12[9] + b2*t12[10];
468 fTnorm[10] = t12[10];
469 fTnorm[11] = 0;
470
471 fTnorm[12] = t12[12] + a2*t12[14];
472 fTnorm[13] = t12[13] + b2*t12[14];
473 fTnorm[14] = t12[14];
474 fTnorm[15] = 1;
475
476 // Scale so that the view volume becomes the canonical one
477 //
478 // Sper = (2/(umax-umin), 2/(vmax-vmin), 1/dproj
479 //
480 Double_t sz = 1./fDproj;
481 Double_t sx = 1./fUVcoord[2];
482 Double_t sy = 1./fUVcoord[3];
483
484 fTnorm[0] *= sx;
485 fTnorm[4] *= sx;
486 fTnorm[8] *= sx;
487 fTnorm[1] *= sy;
488 fTnorm[5] *= sy;
489 fTnorm[9] *= sy;
490 fTnorm[2] *= sz;
491 fTnorm[6] *= sz;
492 fTnorm[10] *= sz;
493 fTnorm[12] *= sx;
494 fTnorm[13] *= sy;
495 fTnorm[14] *= sz;
496}
497
498////////////////////////////////////////////////////////////////////////////////
499/// Define view direction (in spherical coordinates)
500///
501/// Compute transformation matrix from world coordinates
502/// to normalized coordinates (-1 to +1)
503///
504/// Input:
505/// - S(3) - scale factors
506/// - C(3) - centre of scope
507/// - COSPHI - longitude COS
508/// - SINPHI - longitude SIN
509/// - COSTHE - latitude COS (angle between +Z and view direction.)
510/// - SINTHE - latitude SIN
511/// - COSPSI - screen plane rotation angle COS
512/// - SINPSI - screen plane rotation angle SIN
513
515 Double_t cosphi, Double_t sinphi,
516 Double_t costhe, Double_t sinthe,
517 Double_t cospsi, Double_t sinpsi,
518 Double_t *tnorm, Double_t *tback)
519{
520 if (IsPerspective()) {
522 return;
523 }
524 Int_t i, k;
525 Double_t tran[16] /* was [4][4] */, rota[16] /* was [4][4] */;
526 Double_t c1, c2, c3, s1, s2, s3, scalex, scaley, scalez;
527
528 // Parameter adjustments
529 tback -= 5;
530 tnorm -= 5;
531
532 scalex = s[0];
533 scaley = s[1];
534 scalez = s[2];
535
536 //*-*- S E T T R A N S L A T I O N M A T R I X
537 tran[0] = 1 / scalex;
538 tran[1] = 0;
539 tran[2] = 0;
540 tran[3] = -c[0] / scalex;
541
542 tran[4] = 0;
543 tran[5] = 1 / scaley;
544 tran[6] = 0;
545 tran[7] = -c[1] / scaley;
546
547 tran[8] = 0;
548 tran[9] = 0;
549 tran[10] = 1 / scalez;
550 tran[11] = -c[2] / scalez;
551
552 tran[12] = 0;
553 tran[13] = 0;
554 tran[14] = 0;
555 tran[15] = 1;
556
557 //*-*- S E T R O T A T I O N M A T R I X
558 // ( C(PSI) S(PSI) 0) (1 0 0 ) ( C(90+PHI) S(90+PHI) 0)
559 // (-S(PSI) C(PSI) 0) * (0 C(THETA) S(THETA)) * (-S(90+PHI) C(90+PHI) 0)
560 // ( 0 0 1) (0 -S(THETA) C(THETA)) ( 0 0 1)
561 c1 = cospsi;
562 s1 = sinpsi;
563 c2 = costhe;
564 s2 = sinthe;
565 c3 = -sinphi;
566 s3 = cosphi;
567
568 rota[0] = c1*c3 - s1*c2*s3;
569 rota[1] = c1*s3 + s1*c2*c3;
570 rota[2] = s1*s2;
571 rota[3] = 0;
572
573 rota[4] = -s1*c3 - c1* c2*s3;
574 rota[5] = -s1*s3 + c1* c2*c3;
575 rota[6] = c1*s2;
576 rota[7] = 0;
577
578 rota[8] = s2*s3;
579 rota[9] = -s2*c3;
580 rota[10] = c2;
581 rota[11] = 0;
582
583 rota[12] = 0;
584 rota[13] = 0;
585 rota[14] = 0;
586 rota[15] = 1;
587
588 //*-*- F I N D T R A N S F O R M A T I O N M A T R I X
589 for (i = 1; i <= 3; ++i) {
590 for (k = 1; k <= 4; ++k) {
591 tnorm[k + (i << 2)] = rota[(i << 2) - 4]*tran[k - 1] + rota[(i
592 << 2) - 3]*tran[k + 3] + rota[(i << 2) - 2]*tran[k +7]
593 + rota[(i << 2) - 1]*tran[k + 11];
594 }
595 }
596
597 //*-*- S E T B A C K T R A N S L A T I O N M A T R I X
598 tran[0] = scalex;
599 tran[3] = c[0];
600
601 tran[5] = scaley;
602 tran[7] = c[1];
603
604 tran[10] = scalez;
605 tran[11] = c[2];
606
607 //*-*- F I N D B A C K T R A N S F O R M A T I O N
608 for (i = 1; i <= 3; ++i) {
609 for (k = 1; k <= 4; ++k) {
610 tback[k + (i << 2)] = tran[(i << 2) - 4]*rota[(k << 2) - 4] +
611 tran[(i << 2) - 3]*rota[(k << 2) - 3] + tran[(i << 2) -2]
612 *rota[(k << 2) - 2] + tran[(i << 2) - 1]*rota[(k <<2) - 1];
613 }
614 }
615}
616
617////////////////////////////////////////////////////////////////////////////////
618/// Draw the outline of a cube while rotating a 3-d object in the pad.
619
621{
622 TPolyLine3D::DrawOutlineCube(outline,rmin,rmax);
623}
624
625////////////////////////////////////////////////////////////////////////////////
626/// Execute action corresponding to one event.
627
629{
630 ExecuteRotateView(event,px,py);
631}
632
633////////////////////////////////////////////////////////////////////////////////
634/// Execute action corresponding to one event.
635///
636/// This member function is called when a object is clicked with the locator
637///
638/// If Left button clicked in the object area, while the button is kept down
639/// the cube representing the surrounding frame for the corresponding
640/// new latitude and longitude position is drawn.
641
643{
644 static Int_t system, framewasdrawn;
645 static Double_t xrange, yrange, xmin, ymin, longitude1, latitude1, longitude2, latitude2;
646 static Double_t newlatitude, newlongitude, oldlatitude, oldlongitude;
647 Double_t dlatitude, dlongitude, x, y;
648 Int_t irep = 0;
649 Double_t psideg;
650 Bool_t opaque = gPad->OpaqueMoving();
651
652 // All coordinates transformation are from absolute to relative
653 if (!gPad->IsEditable()) return;
654 gPad->AbsCoordinates(kTRUE);
655
656 switch (event) {
657
658 case kKeyPress :
659 fChanged = kTRUE;
660 MoveViewCommand(Char_t(px), py);
661 break;
662
663 case kMouseMotion:
664 gPad->SetCursor(kRotate);
665 break;
666
667 case kButton1Down:
668
669 // remember position of the cube
670 xmin = gPad->GetX1();
671 ymin = gPad->GetY1();
672 xrange = gPad->GetX2() - xmin;
673 yrange = gPad->GetY2() - ymin;
674 x = gPad->PixeltoX(px);
675 y = gPad->PixeltoY(py);
676 system = GetSystem();
677 framewasdrawn = 0;
678 if (system == kCARTESIAN || system == kPOLAR || system == 11 || IsPerspective()) {
679 longitude1 = 180*(x-xmin)/xrange;
680 latitude1 = 90*(y-ymin)/yrange;
681 } else {
682 latitude1 = 90*(x-xmin)/xrange;
683 longitude1 = 180*(y-ymin)/yrange;
684 }
685 newlongitude = oldlongitude = -90 - gPad->GetPhi();
686 newlatitude = oldlatitude = 90 - gPad->GetTheta();
687 psideg = GetPsi();
688
689 // if outline isn't set, make it look like a cube
690 if(!fOutline)
692 break;
693
694 case kButton1Motion:
695 {
696 // draw the surrounding frame for the current mouse position
697 // first: Erase old frame
698 fChanged = kTRUE;
699 if (framewasdrawn && !opaque) fOutline->Paint();
700 framewasdrawn = 1;
701 x = gPad->PixeltoX(px);
702 y = gPad->PixeltoY(py);
703 if (system == kCARTESIAN || system == kPOLAR || system == 11 || IsPerspective()) {
704 longitude2 = 180*(x-xmin)/xrange;
705 latitude2 = 90*(y-ymin)/yrange;
706 } else {
707 latitude2 = 90*(x-xmin)/xrange;
708 longitude2 = 180*(y-ymin)/yrange;
709 }
710 dlongitude = longitude2 - longitude1;
711 dlatitude = latitude2 - latitude1;
712 newlatitude = oldlatitude + dlatitude;
713 newlongitude = oldlongitude - dlongitude;
714 psideg = GetPsi();
715 ResetView(newlongitude, newlatitude, psideg, irep);
716 if (!opaque) {
717 fOutline->Paint();
718 } else {
719 psideg = GetPsi();
720 SetView(newlongitude, newlatitude, psideg, irep);
721 gPad->SetPhi(-90-newlongitude);
722 gPad->SetTheta(90-newlatitude);
723 gPad->Modified(kTRUE);
724 }
725
726 break;
727 }
728 case kButton1Up:
729 if (gROOT->IsEscaped()) {
730 gROOT->SetEscape(kFALSE);
731 if (opaque) {
732 psideg = GetPsi();
733 SetView(oldlongitude, oldlatitude, psideg, irep);
734 gPad->SetPhi(-90-oldlongitude);
735 gPad->SetTheta(90-oldlatitude);
736 gPad->Modified(kTRUE);
737 }
738 break;
739 }
740
741 // Temporary fix for 2D drawing problems on pad. fOutline contains
742 // a TPolyLine3D object for the rotation box. This will be painted
743 // through a newly created TView3Der3DPad instance, which is left
744 // behind on pad. This remaining creates 2D drawing problems.
745 //
746 // This is a TEMPORARY fix - will be removed when proper multiple viewers
747 // on pad problems are resolved.
748 if (gPad) {
749 TVirtualViewer3D *viewer = gPad->GetViewer3D();
750 if (viewer && !strcmp(viewer->IsA()->GetName(),"TView3Der3DPad")) {
751 gPad->ReleaseViewer3D();
752 delete viewer;
753 }
754 }
755 // End fix
756
757 // Recompute new view matrix and redraw
758 psideg = GetPsi();
759 SetView(newlongitude, newlatitude, psideg, irep);
760 gPad->SetPhi(-90-newlongitude);
761 gPad->SetTheta(90-newlatitude);
762 gPad->Modified(kTRUE);
763
764 // Set line color, style and width
765 gVirtualX->SetLineColor(-1);
766 gVirtualX->SetLineStyle(-1);
767 gVirtualX->SetLineWidth(-1);
768 break;
769 }
770
771 // set back to default transformation mode
772 gPad->AbsCoordinates(kFALSE);
773}
774
775////////////////////////////////////////////////////////////////////////////////
776/// Find Z component of NORMAL in normalized coordinates.
777///
778/// Input:
779/// - X - X-component of NORMAL
780/// - Y - Y-component of NORMAL
781/// - Z - Z-component of NORMAL
782///
783/// Output:
784/// - ZN - Z-component of NORMAL in normalized coordinates
785
787{
788 zn = x*(fTN[1] * fTN[6] - fTN[2] * fTN[5]) + y*(fTN[2] * fTN[4] -
789 fTN[0] * fTN[6]) + z*(fTN[0] * fTN[5] - fTN[1] * fTN[4]);
790}
791
792////////////////////////////////////////////////////////////////////////////////
793/// Find critical PHI sectors.
794///
795/// Input:
796/// - IOPT - options:
797/// 1. from BACK to FRONT 'BF'
798/// 2. from FRONT to BACK 'FB'
799/// - KPHI - number of phi sectors
800/// - APHI(*) - PHI separators (modified internally)
801///
802/// Output:
803/// - IPHI1 - initial sector
804/// - IPHI2 - final sector
805
806void TView3D::FindPhiSectors(Int_t iopt, Int_t &kphi, Double_t *aphi, Int_t &iphi1, Int_t &iphi2)
807{
808 Int_t iphi[2], i, k;
809 Double_t dphi;
810 Double_t x1, x2, z1, z2, phi1, phi2;
811
812 // Parameter adjustments
813 --aphi;
814
815 if (aphi[kphi + 1] == aphi[1]) aphi[kphi + 1] += 360;
816 dphi = TMath::Abs(aphi[kphi + 1] - aphi[1]);
817 if (dphi != 360) {
818 aphi[kphi + 2] = (aphi[1] + aphi[kphi + 1]) / (float)2. + 180;
819 aphi[kphi + 3] = aphi[1] + 360;
820 kphi += 2;
821 }
822
823 //*-*- F I N D C R I T I C A L S E C T O R S
824 k = 0;
825 for (i = 1; i <= kphi; ++i) {
826 phi1 = kRad*aphi[i];
827 phi2 = kRad*aphi[i + 1];
828 x1 = fTN[0]*TMath::Cos(phi1) + fTN[1]*TMath::Sin(phi1);
829 x2 = fTN[0]*TMath::Cos(phi2) + fTN[1]*TMath::Sin(phi2);
830 if (x1 >= 0 && x2 > 0) continue;
831 if (x1 <= 0 && x2 < 0) continue;
832 ++k;
833 if (k == 3) break;
834 iphi[k - 1] = i;
835 }
836 if (k != 2) {
837 Error("FindPhiSectors", "something strange: num. of critical sector not equal 2");
838 iphi1 = 1;
839 iphi2 = 2;
840 return;
841 }
842
843 //*-*- F I N D O R D E R O F C R I T I C A L S E C T O R S
844 phi1 = kRad*(aphi[iphi[0]] + aphi[iphi[0] + 1]) / (float)2.;
845 phi2 = kRad*(aphi[iphi[1]] + aphi[iphi[1] + 1]) / (float)2.;
846 z1 = fTN[8]*TMath::Cos(phi1) + fTN[9]*TMath::Sin(phi1);
847 z2 = fTN[8]*TMath::Cos(phi2) + fTN[9]*TMath::Sin(phi2);
848 if ((z1 <= z2 && iopt == 1) || (z1 > z2 && iopt == 2)) {
849 iphi1 = iphi[0];
850 iphi2 = iphi[1];
851 } else {
852 iphi1 = iphi[1];
853 iphi2 = iphi[0];
854 }
855}
856
857////////////////////////////////////////////////////////////////////////////////
858/// Find critical THETA sectors for given PHI sector.
859///
860/// Input:
861/// - IOPT - options:
862/// 1. from BACK to FRONT 'BF'
863/// 2. from FRONT to BACK 'FB'
864/// - PHI - PHI sector
865/// - KTH - number of THETA sectors
866/// - ATH(*) - THETA separators (modified internally)
867///
868/// Output:
869/// - ITH1 - initial sector
870/// - ITH2 - final sector
871
872void TView3D::FindThetaSectors(Int_t iopt, Double_t phi, Int_t &kth, Double_t *ath, Int_t &ith1, Int_t &ith2)
873{
874 Int_t i, k, ith[2];
875 Double_t z1, z2, cosphi, sinphi, tncons, th1, th2, dth;
876
877 // Parameter adjustments
878 --ath;
879
880 // Function Body
881 dth = TMath::Abs(ath[kth + 1] - ath[1]);
882 if (dth != 360) {
883 ath[kth + 2] = 0.5*(ath[1] + ath[kth + 1]) + 180;
884 ath[kth + 3] = ath[1] + 360;
885 kth += 2;
886 }
887
888 //*-*- F I N D C R I T I C A L S E C T O R S
889 cosphi = TMath::Cos(phi*kRad);
890 sinphi = TMath::Sin(phi*kRad);
891 k = 0;
892 for (i = 1; i <= kth; ++i) {
893 th1 = kRad*ath[i];
894 th2 = kRad*ath[i + 1];
895 FindNormal(TMath::Cos(th1)*cosphi, TMath::Cos(th1)*sinphi, -TMath::Sin(th1), z1);
896 FindNormal(TMath::Cos(th2)*cosphi, TMath::Cos(th2)*sinphi, -TMath::Sin(th2), z2);
897 if (z1 >= 0 && z2 > 0) continue;
898 if (z1 <= 0 && z2 < 0) continue;
899 ++k;
900 if (k == 3) break;
901 ith[k - 1] = i;
902 }
903 if (k != 2) {
904 Error("FindThetaSectors", "Something strange: num. of critical sectors not equal 2");
905 ith1 = 1;
906 ith2 = 2;
907 return;
908 }
909
910 //*-*- F I N D O R D E R O F C R I T I C A L S E C T O R S
911 tncons = fTN[8]*TMath::Cos(phi*kRad) + fTN[9]*TMath::Sin(phi*kRad);
912 th1 = kRad*(ath[ith[0]] + ath[ith[0] + 1]) / (float)2.;
913 th2 = kRad*(ath[ith[1]] + ath[ith[1] + 1]) / (float)2.;
914 z1 = tncons*TMath::Sin(th1) + fTN[10]*TMath::Cos(th1);
915 z2 = tncons*TMath::Sin(th2) + fTN[10]*TMath::Cos(th2);
916 if ((z1 <= z2 && iopt == 1) || (z1 > z2 && iopt == 2)) {
917 ith1 = ith[0];
918 ith2 = ith[1];
919 } else {
920 ith1 = ith[1];
921 ith2 = ith[0];
922 }
923}
924
925////////////////////////////////////////////////////////////////////////////////
926/// Find centre of a MIN-MAX scope and scale factors
927///
928/// Output:
929/// - SCALE(3) - scale factors
930/// - CENTER(3) - centre
931/// - IREP - reply (-1 if error in min-max)
932
933void TView3D::FindScope(Double_t *scale, Double_t *center, Int_t &irep)
934{
935 irep = 0;
936 Double_t sqrt3 = 0.5*TMath::Sqrt(3.0);
937
938 for (Int_t i = 0; i < 3; i++) {
939 if (fRmin[i] >= fRmax[i]) { irep = -1; return;}
940 scale[i] = sqrt3*(fRmax[i] - fRmin[i]);
941 center[i] = 0.5*(fRmax[i] + fRmin[i]);
942 }
943}
944
945////////////////////////////////////////////////////////////////////////////////
946/// Return distance to axis from point px,py.
947///
948/// Algorithm:
949///
950/// ~~~ {.cpp}
951/// A(x1,y1) P B(x2,y2)
952/// ------------------------------------------------
953/// I
954/// I
955/// I
956/// I
957/// M(x,y)
958///
959/// Let us call a = distance AM A=a**2
960/// b = distance BM B=b**2
961/// c = distance AB C=c**2
962/// d = distance PM D=d**2
963/// u = distance AP U=u**2
964/// v = distance BP V=v**2 c = u + v
965///
966/// D = A - U
967/// D = B - V = B -(c-u)**2
968/// ==> u = (A -B +C)/2c
969/// ~~~
970
972{
973 Double_t x1,y1,x2,y2;
974 Double_t x = px;
975 Double_t y = py;
976 ratio = 0;
977
978 if (fSystem != kCARTESIAN) return 9998; // only implemented for Cartesian coordinates
979 if (axis == 1) {
980 x1 = gPad->XtoAbsPixel(fX1[0]);
981 y1 = gPad->YtoAbsPixel(fX1[1]);
982 x2 = gPad->XtoAbsPixel(fX2[0]);
983 y2 = gPad->YtoAbsPixel(fX2[1]);
984 } else if (axis == 2) {
985 x1 = gPad->XtoAbsPixel(fY1[0]);
986 y1 = gPad->YtoAbsPixel(fY1[1]);
987 x2 = gPad->XtoAbsPixel(fY2[0]);
988 y2 = gPad->YtoAbsPixel(fY2[1]);
989 } else {
990 x1 = gPad->XtoAbsPixel(fZ1[0]);
991 y1 = gPad->YtoAbsPixel(fZ1[1]);
992 x2 = gPad->XtoAbsPixel(fZ2[0]);
993 y2 = gPad->YtoAbsPixel(fZ2[1]);
994 }
995 Double_t xx1 = x - x1;
996 Double_t xx2 = x - x2;
997 Double_t x1x2 = x1 - x2;
998 Double_t yy1 = y - y1;
999 Double_t yy2 = y - y2;
1000 Double_t y1y2 = y1 - y2;
1001 Double_t a = xx1*xx1 + yy1*yy1;
1002 Double_t b = xx2*xx2 + yy2*yy2;
1003 Double_t c = x1x2*x1x2 + y1y2*y1y2;
1004 if (c <= 0) return 9999;
1006 Double_t u = (a - b + c)/(2*v);
1007 Double_t d = TMath::Abs(a - u*u);
1008
1009 Int_t dist = Int_t(TMath::Sqrt(d) - 0.5);
1010 ratio = u/v;
1011 return dist;
1012}
1013
1014////////////////////////////////////////////////////////////////////////////////
1015/// Get maximum view extent.
1016
1018{
1019 Double_t dx = 0.5*(fRmax[0]-fRmin[0]);
1020 Double_t dy = 0.5*(fRmax[1]-fRmin[1]);
1021 Double_t dz = 0.5*(fRmax[2]-fRmin[2]);
1022 Double_t extent = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
1023 return extent;
1024}
1025
1026////////////////////////////////////////////////////////////////////////////////
1027/// Get Range function.
1028
1030{
1031 for (Int_t i = 0; i < 3; max[i] = fRmax[i], min[i] = fRmin[i], i++) { }
1032}
1033
1034////////////////////////////////////////////////////////////////////////////////
1035/// Get Range function.
1036
1038{
1039 for (Int_t i = 0; i < 3; max[i] = fRmax[i], min[i] = fRmin[i], i++) { }
1040}
1041
1042////////////////////////////////////////////////////////////////////////////////
1043/// Get current window extent.
1044
1046{
1047 u0 = fUVcoord[0];
1048 v0 = fUVcoord[1];
1049 du = fUVcoord[2];
1050 dv = fUVcoord[3];
1051}
1052
1053////////////////////////////////////////////////////////////////////////////////
1054/// Check if point is clipped in perspective view.
1055
1057{
1058 if (TMath::Abs(p[0])>p[2]) return kTRUE;
1059 if (TMath::Abs(p[1])>p[2]) return kTRUE;
1060 return kFALSE;
1061}
1062
1063////////////////////////////////////////////////////////////////////////////////
1064/// Transfer point from normalized to world coordinates.
1065///
1066/// Input:
1067/// - PN(3) - point in world coordinate system
1068/// - PW(3) - point in normalized coordinate system
1069
1071{
1072 Float_t x = pn[0], y = pn[1], z = pn[2];
1073 pw[0] = fTback[0]*x + fTback[1]*y + fTback[2]*z + fTback[3];
1074 pw[1] = fTback[4]*x + fTback[5]*y + fTback[6]*z + fTback[7];
1075 pw[2] = fTback[8]*x + fTback[9]*y + fTback[10]*z + fTback[11];
1076}
1077
1078////////////////////////////////////////////////////////////////////////////////
1079/// Transfer point from normalized to world coordinates.
1080///
1081/// Input:
1082/// - PN(3) - point in world coordinate system
1083/// - PW(3) - point in normalized coordinate system
1084
1086{
1087 Double_t x = pn[0], y = pn[1], z = pn[2];
1088 pw[0] = fTback[0]*x + fTback[1]*y + fTback[2]*z + fTback[3];
1089 pw[1] = fTback[4]*x + fTback[5]*y + fTback[6]*z + fTback[7];
1090 pw[2] = fTback[8]*x + fTback[9]*y + fTback[10]*z + fTback[11];
1091}
1092
1093////////////////////////////////////////////////////////////////////////////////
1094/// Transfer vector of NORMAL from word to normalized coordinates.
1095///
1096/// Input:
1097/// - PW(3) - vector of NORMAL in word coordinate system
1098/// - PN(3) - vector of NORMAL in normalized coordinate system
1099
1101{
1102 Double_t x, y, z, a1, a2, a3, b1, b2, b3, c1, c2, c3;
1103
1104 x = pw[0];
1105 y = pw[1];
1106 z = pw[2];
1107 a1 = fTnorm[0];
1108 a2 = fTnorm[1];
1109 a3 = fTnorm[2];
1110 b1 = fTnorm[4];
1111 b2 = fTnorm[5];
1112 b3 = fTnorm[6];
1113 c1 = fTnorm[8];
1114 c2 = fTnorm[9];
1115 c3 = fTnorm[10];
1116 pn[0] = x*(b2*c3 - b3*c2) + y*(b3*c1 - b1*c3) + z*(b1*c2 - b2*c1);
1117 pn[1] = x*(c2*a3 - c3*a2) + y*(c3*a1 - c1*a3) + z*(c1*a2 - c2*a1);
1118 pn[2] = x*(a2*b3 - a3*b2) + y*(a3*b1 - a1*b3) + z*(a1*b2 - a2*b1);
1119}
1120
1121////////////////////////////////////////////////////////////////////////////////
1122/// Transfer vector of NORMAL from word to normalized coordinates.
1123///
1124/// Input:
1125/// - PW(3) - vector of NORMAL in word coordinate system
1126/// - PN(3) - vector of NORMAL in normalized coordinate system
1127
1129{
1130 Double_t x, y, z, a1, a2, a3, b1, b2, b3, c1, c2, c3;
1131
1132 x = pw[0];
1133 y = pw[1];
1134 z = pw[2];
1135 a1 = fTnorm[0];
1136 a2 = fTnorm[1];
1137 a3 = fTnorm[2];
1138 b1 = fTnorm[4];
1139 b2 = fTnorm[5];
1140 b3 = fTnorm[6];
1141 c1 = fTnorm[8];
1142 c2 = fTnorm[9];
1143 c3 = fTnorm[10];
1144 pn[0] = x*(b2*c3 - b3*c2) + y*(b3*c1 - b1*c3) + z*(b1*c2 - b2*c1);
1145 pn[1] = x*(c2*a3 - c3*a2) + y*(c3*a1 - c1*a3) + z*(c1*a2 - c2*a1);
1146 pn[2] = x*(a2*b3 - a3*b2) + y*(a3*b1 - a1*b3) + z*(a1*b2 - a2*b1);
1147}
1148
1149////////////////////////////////////////////////////////////////////////////////
1150/// Set the correct window size for lego and surface plots.
1151///
1152/// Set the correct window size for lego and surface plots.
1153/// And draw the background if necessary.
1154///
1155/// Input parameters:
1156/// - RBACK : Background colour
1157
1159{
1160 Int_t i, k;
1161 Double_t x, y, z, r1, r2, r3, xx, yy, smax[2];
1162 Double_t xgraf[6], ygraf[6];
1163
1164 for (i = 1; i <= 2; ++i) {
1165 smax[i - 1] = fTnorm[(i << 2) - 1];
1166 for (k = 1; k <= 3; ++k) {
1167 if (fTnorm[k + (i << 2) - 5] < 0) {
1168 smax[i - 1] += fTnorm[k + (i << 2) - 5]*fRmin[k-1];
1169 } else {
1170 smax[i - 1] += fTnorm[k + (i << 2) - 5]*fRmax[k-1];
1171 }
1172 }
1173 }
1174
1175 //*-*- Compute x,y range
1176 Double_t xmin = -smax[0];
1177 Double_t xmax = smax[0];
1178 Double_t ymin = -smax[1];
1179 Double_t ymax = smax[1];
1180 Double_t dx = xmax-xmin;
1181 Double_t dy = ymax-ymin;
1182 Double_t dxr = dx/(1 - gPad->GetLeftMargin() - gPad->GetRightMargin());
1183 Double_t dyr = dy/(1 - gPad->GetBottomMargin() - gPad->GetTopMargin());
1184
1185 // Range() could change the size of the pad pixmap and therefore should
1186 // be called before the other paint routines
1187 gPad->Range(xmin - dxr*gPad->GetLeftMargin(),
1188 ymin - dyr*gPad->GetBottomMargin(),
1189 xmax + dxr*gPad->GetRightMargin(),
1190 ymax + dyr*gPad->GetTopMargin());
1191 gPad->RangeAxis(xmin, ymin, xmax, ymax);
1192
1193 //*-*- Draw the background if necessary
1194 if (rback > 0) {
1195 r1 = -1;
1196 r2 = -1;
1197 r3 = -1;
1198 xgraf[0] = -smax[0];
1199 xgraf[1] = -smax[0];
1200 xgraf[2] = -smax[0];
1201 xgraf[3] = -smax[0];
1202 xgraf[4] = smax[0];
1203 xgraf[5] = smax[0];
1204 ygraf[0] = -smax[1];
1205 ygraf[1] = smax[1];
1206 ygraf[2] = -smax[1];
1207 ygraf[3] = smax[1];
1208 ygraf[5] = smax[1];
1209 ygraf[4] = -smax[1];
1210 for (i = 1; i <= 8; ++i) {
1211 x = 0.5*((1 - r1)*fRmin[0] + (r1 + 1)*fRmax[0]);
1212 y = 0.5*((1 - r2)*fRmin[1] + (r2 + 1)*fRmax[1]);
1213 z = 0.5*((1 - r3)*fRmin[2] + (r3 + 1)*fRmax[2]);
1214 xx = fTnorm[0]*x + fTnorm[1]*y + fTnorm[2]*z + fTnorm[3];
1215 yy = fTnorm[4]*x + fTnorm[5]*y + fTnorm[6]*z + fTnorm[7];
1216 if (TMath::Abs(xx - xgraf[1]) <= 1e-4) {
1217 if (ygraf[1] >= yy) ygraf[1] = yy;
1218 if (ygraf[2] <= yy) ygraf[2] = yy;
1219 }
1220 if (TMath::Abs(xx - xgraf[5]) <= 1e-4) {
1221 if (ygraf[5] >= yy) ygraf[5] = yy;
1222 if (ygraf[4] <= yy) ygraf[4] = yy;
1223 }
1224 if (TMath::Abs(yy - ygraf[0]) <= 1e-4) xgraf[0] = xx;
1225 if (TMath::Abs(yy - ygraf[3]) <= 1e-4) xgraf[3] = xx;
1226 r1 = -r1;
1227 if (i % 2 == 0) r2 = -r2;
1228 if (i >= 4) r3 = 1;
1229 }
1230 gPad->PaintFillArea(6, xgraf, ygraf);
1231 }
1232}
1233
1234////////////////////////////////////////////////////////////////////////////////
1235/// Store axis coordinates in the NDC system.
1236
1237void TView3D::SetAxisNDC(const Double_t *x1, const Double_t *x2, const Double_t *y1, const Double_t *y2, const Double_t *z1, const Double_t *z2)
1238{
1239 for (Int_t i=0;i<3;i++) {
1240 fX1[i] = x1[i];
1241 fX2[i] = x2[i];
1242 fY1[i] = y1[i];
1243 fY2[i] = y2[i];
1244 fZ1[i] = z1[i];
1245 fZ2[i] = z2[i];
1246 }
1247}
1248
1249////////////////////////////////////////////////////////////////////////////////
1250/// Set default viewing window.
1251
1253{
1254 if (!gPad) return;
1255 Double_t screen_factor = 1.;
1256 Double_t du, dv;
1257 Double_t extent = GetExtent();
1258 fDview = 3*extent;
1259 fDproj = 0.5*extent;
1260
1261 // width in pixels
1262 fUpix = gPad->GetWw()*gPad->GetAbsWNDC();
1263
1264 // height in pixels
1265 fVpix = gPad->GetWh()*gPad->GetAbsHNDC();
1266 du = 0.5*screen_factor*fDproj;
1267 dv = du*fVpix/fUpix; // keep aspect ratio
1268 SetWindow(0, 0, du, dv);
1269}
1270
1271////////////////////////////////////////////////////////////////////////////////
1272/// This is a function which creates default outline.
1273///
1274/// ~~~ {.cpp}
1275/// x = fRmin[0] X = fRmax[0]
1276/// y = fRmin[1] Y = fRmax[1]
1277/// z = fRmin[2] Z = fRmax[2]
1278///
1279/// (x,Y,Z) +---------+ (X,Y,Z)
1280/// / /|
1281/// / / |
1282/// / / |
1283/// (x,y,Z) +---------+ |
1284/// | | + (X,Y,z)
1285/// | | /
1286/// | | /
1287/// | |/
1288/// +---------+
1289/// (x,y,z) (X,y,z)
1290/// ~~~
1291
1293{
1294 if (!fOutline) {
1296 fOutline = new TList();
1297 }
1299}
1300
1301////////////////////////////////////////////////////////////////////////////////
1302/// Set the parallel option (default).
1303
1305{
1306 if (!IsPerspective()) return;
1308 Int_t irep;
1310}
1311
1312////////////////////////////////////////////////////////////////////////////////
1313/// Set perspective option.
1314
1316{
1317 if (IsPerspective()) return;
1319 Int_t irep;
1322}
1323
1324////////////////////////////////////////////////////////////////////////////////
1325/// Set Range function.
1326
1327void TView3D::SetRange(const Double_t *min, const Double_t *max)
1328{
1329 Int_t irep;
1330 for (Int_t i = 0; i < 3; fRmax[i] = max[i], fRmin[i] = min[i], i++) { }
1333 if(irep < 0)
1334 Error("SetRange", "problem setting view");
1336}
1337
1338////////////////////////////////////////////////////////////////////////////////
1339/// Set 3-D View range.
1340///
1341/// Input:
1342/// - x0, y0, z0 are minimum coordinates
1343/// - x1, y1, z1 are maximum coordinates
1344///
1345/// - flag values are:
1346/// - 0 (set always) <- default
1347/// - 1 (shrink view)
1348/// - 2 (expand view)
1349
1351{
1352 Double_t rmax[3], rmin[3];
1353
1354 switch (flag) {
1355 case 2: // expand view
1356 GetRange(rmin, rmax);
1357 rmin[0] = x0 < rmin[0] ? x0 : rmin[0];
1358 rmin[1] = y0 < rmin[1] ? y0 : rmin[1];
1359 rmin[2] = z0 < rmin[2] ? z0 : rmin[2];
1360 rmax[0] = x1 > rmax[0] ? x1 : rmax[0];
1361 rmax[1] = y1 > rmax[1] ? y1 : rmax[1];
1362 rmax[2] = z1 > rmax[2] ? z1 : rmax[2];
1363 break;
1364
1365 case 1: // shrink view
1366 GetRange(rmin, rmax);
1367 rmin[0] = x0 > rmin[0] ? x0 : rmin[0];
1368 rmin[1] = y0 > rmin[1] ? y0 : rmin[1];
1369 rmin[2] = z0 > rmin[2] ? z0 : rmin[2];
1370 rmax[0] = x1 < rmax[0] ? x1 : rmax[0];
1371 rmax[1] = y1 < rmax[1] ? y1 : rmax[1];
1372 rmax[2] = z1 < rmax[2] ? z1 : rmax[2];
1373 break;
1374
1375 default:
1376 rmin[0] = x0; rmax[0] = x1;
1377 rmin[1] = y0; rmax[1] = y1;
1378 rmin[2] = z0; rmax[2] = z1;
1379 }
1380 SetRange(rmin, rmax);
1381}
1382
1383////////////////////////////////////////////////////////////////////////////////
1384/// Set viewing window.
1385
1387{
1388 fUVcoord[0] = u0;
1389 fUVcoord[1] = v0;
1390 fUVcoord[2] = du;
1391 fUVcoord[3] = dv;
1392}
1393
1394////////////////////////////////////////////////////////////////////////////////
1395/// Set view parameters.
1396
1397void TView3D::SetView(Double_t longitude, Double_t latitude, Double_t psi, Int_t &irep)
1398{
1399 ResetView(longitude, latitude, psi, irep);
1400}
1401
1402////////////////////////////////////////////////////////////////////////////////
1403/// Recompute window for perspective view.
1404
1406{
1407 if (!IsPerspective()) return;
1408 Double_t upix = fUpix;
1409 Double_t vpix = fVpix;
1410
1411 // width in pixels
1412 fUpix = gPad->GetWw()*gPad->GetAbsWNDC();
1413
1414 // height in pixels
1415 fVpix = gPad->GetWh()*gPad->GetAbsHNDC();
1416 Double_t u0 = fUVcoord[0]*fUpix/upix;
1417 Double_t v0 = fUVcoord[1]*fVpix/vpix;
1418 Double_t du = fUVcoord[2]*fUpix/upix;
1419 Double_t dv = fUVcoord[3]*fVpix/vpix;
1420 SetWindow(u0, v0, du, dv);
1422}
1423
1424////////////////////////////////////////////////////////////////////////////////
1425/// Set view direction (in spherical coordinates).
1426///
1427/// Input
1428/// - PHI - longitude
1429/// - THETA - latitude (angle between +Z and view direction)
1430/// - PSI - rotation in screen plane
1431///
1432/// Output:
1433/// - IREP - reply (-1 if error in min-max)
1434///
1435/// Errors: error in min-max scope
1436
1437void TView3D::ResetView(Double_t longitude, Double_t latitude, Double_t psi, Int_t &irep)
1438{
1439 Double_t scale[3], centre[3];
1440 Double_t c1, c2, c3, s1, s2, s3;
1441
1442 //*-*- F I N D C E N T E R O F S C O P E A N D
1443 //*-*- S C A L E F A C T O R S
1444 FindScope(scale, centre, irep);
1445 if (irep < 0) {
1446 Error("ResetView", "Error in min-max scope");
1447 return;
1448 }
1449
1450 //*-*- S E T T R A N S F O R M A T I O N M A T R I C E S
1451 fLongitude = longitude;
1452 fPsi = psi;
1453 fLatitude = latitude;
1454
1455 if (IsPerspective()) {
1457 return;
1458 }
1459
1460 c1 = TMath::Cos(longitude*kRad);
1461 s1 = TMath::Sin(longitude*kRad);
1462 c2 = TMath::Cos(latitude*kRad);
1463 s2 = TMath::Sin(latitude*kRad);
1464 c3 = TMath::Cos(psi*kRad);
1465 s3 = TMath::Sin(psi*kRad);
1466 DefineViewDirection(scale, centre, c1, s1, c2, s2, c3, s3, fTnorm, fTback);
1467 c3 = 1;
1468 s3 = 0;
1469 DefineViewDirection(scale, centre, c1, s1, c2, s2, c3, s3, fTN, fTB);
1470}
1471
1472////////////////////////////////////////////////////////////////////////////////
1473/// Transfer point from world to normalized coordinates.
1474///
1475/// Input:
1476/// - PW(3) - point in world coordinate system
1477/// - PN(3) - point in normalized coordinate system
1478
1480{
1481 Float_t x = pw[0], y = pw[1], z = pw[2];
1482
1483 // perspective view
1484 if (IsPerspective()) {
1485 for (Int_t i=0; i<3; i++) {
1486 pn[i] = fTnorm[i]*x + fTnorm[i+4]*y + fTnorm[i+8]*z + fTnorm[i+12];
1487 }
1488 if (pn[2]>0) {
1489 pn[0] /= pn[2];
1490 pn[1] /= pn[2];
1491 } else {
1492 pn[0] *= 1000.;
1493 pn[1] *= 1000.;
1494 }
1495 return;
1496 }
1497
1498 // parallel view
1499 pn[0] = fTnorm[0]*x + fTnorm[1]*y + fTnorm[2]*z + fTnorm[3];
1500 pn[1] = fTnorm[4]*x + fTnorm[5]*y + fTnorm[6]*z + fTnorm[7];
1501 pn[2] = fTnorm[8]*x + fTnorm[9]*y + fTnorm[10]*z + fTnorm[11];
1502}
1503
1504////////////////////////////////////////////////////////////////////////////////
1505/// Transfer point from world to normalized coordinates.
1506///
1507/// Input:
1508/// - PW(3) - point in world coordinate system
1509/// - PN(3) - point in normalized coordinate system
1510
1512{
1513 Double_t x = pw[0], y = pw[1], z = pw[2];
1514
1515 // perspective view
1516 if (IsPerspective()) {
1517 for (Int_t i=0; i<3; i++) {
1518 pn[i] = fTnorm[i]*x + fTnorm[i+4]*y + fTnorm[i+8]*z + fTnorm[i+12];
1519 }
1520 if (pn[2]>0) {
1521 pn[0] /= pn[2];
1522 pn[1] /= pn[2];
1523 } else {
1524 pn[0] *= 1000.;
1525 pn[1] *= 1000.;
1526 }
1527 return;
1528 }
1529
1530 // parallel view
1531 pn[0] = fTnorm[0]*x + fTnorm[1]*y + fTnorm[2]*z + fTnorm[3];
1532 pn[1] = fTnorm[4]*x + fTnorm[5]*y + fTnorm[6]*z + fTnorm[7];
1533 pn[2] = fTnorm[8]*x + fTnorm[9]*y + fTnorm[10]*z + fTnorm[11];
1534}
1535
1536////////////////////////////////////////////////////////////////////////////////
1537/// Force the current pad to be updated.
1538
1540{
1541 TVirtualPad *thisPad = pad;
1542 if (!thisPad) thisPad = gPad;
1543 if (thisPad) {
1544#ifdef R__HAS_COCOA
1545 thisPad->AbsCoordinates(kFALSE);
1546#endif
1547 thisPad->Modified();
1548 thisPad->Update();
1549 }
1550}
1551
1552////////////////////////////////////////////////////////////////////////////////
1553/// API to rotate view and adjust the pad provided it the current one.
1554
1556{
1557 Int_t iret;
1558 Double_t p = phi;
1559 Double_t t = theta;
1560 SetView(p, t, 0, iret);
1561
1562 // Adjust current pad too
1563 TVirtualPad *thisPad = pad;
1564 if (!thisPad) thisPad = gPad;
1565 if (thisPad) {
1566 thisPad->SetPhi(-90-p);
1567 thisPad->SetTheta(90-t);
1568 thisPad->Modified();
1569 thisPad->Update();
1570 }
1571}
1572
1573////////////////////////////////////////////////////////////////////////////////
1574/// Set to side view.
1575
1577{
1578 RotateView(0,90.0,pad);
1579}
1580
1581////////////////////////////////////////////////////////////////////////////////
1582/// Set to front view.
1583
1585{
1586 RotateView(270.0,90.0,pad);
1587}
1588
1589////////////////////////////////////////////////////////////////////////////////
1590/// Set to top view.
1591
1593{
1594 RotateView(270.0,0.0,pad);
1595}
1596
1597////////////////////////////////////////////////////////////////////////////////
1598/// Turn on /off 3D axis.
1599
1601{
1603}
1604
1605////////////////////////////////////////////////////////////////////////////////
1606/// Turn on /off the interactive option to
1607/// Zoom / Move / Change attributes of 3D axis correspond this view.
1608
1610{
1612}
1613
1614////////////////////////////////////////////////////////////////////////////////
1615/// Adjust all sides of view in respect of the biggest one.
1616
1618{
1619 Double_t min[3],max[3];
1620 GetRange(min,max);
1621 int i;
1622 Double_t maxSide = 0;
1623 // Find the largest side
1624 for (i=0;i<3; i++) maxSide = TMath::Max(maxSide,max[i]-min[i]);
1625 //Adjust scales:
1626 for (i=0;i<3; i++) max[i] += maxSide - (max[i]-min[i]);
1627 SetRange(min,max);
1628
1629 AdjustPad(pad);
1630}
1631
1632////////////////////////////////////////////////////////////////////////////////
1633/// Move view into the center of the scene.
1634
1636{
1637 Double_t min[3],max[3];
1638 GetRange(min,max);
1639 int i;
1640 for (i=0;i<3; i++) {
1641 if (max[i] > 0) min[i] = -max[i];
1642 else max[i] = -min[i];
1643 }
1644 SetRange(min,max);
1645 AdjustPad(pad);
1646}
1647
1648////////////////////////////////////////////////////////////////////////////////
1649/// unZOOM this view.
1650
1652{
1653 if (TMath::Abs(unZoomFactor) < 0.001) return;
1654 ZoomView(pad,1./unZoomFactor);
1655}
1656
1657////////////////////////////////////////////////////////////////////////////////
1658/// ZOOM this view.
1659
1661{
1662 if (TMath::Abs(zoomFactor) < 0.001) return;
1663 Double_t min[3],max[3];
1664 GetRange(min,max);
1665 int i;
1666 for (i=0;i<3; i++) {
1667 // Find center
1668 Double_t c = (max[i]+min[i])/2;
1669 // Find a new size
1670 Double_t s = (max[i]-min[i])/(2*zoomFactor);
1671 // Set a new size
1672 max[i] = c + s;
1673 min[i] = c - s;
1674 }
1675 SetRange(min,max);
1676 AdjustPad(pad);
1677}
1678
1679////////////////////////////////////////////////////////////////////////////////
1680/// Move focus to a different box position and extent in nsteps. Perform
1681/// rotation with dlat,dlong,dpsi at each step.
1682
1684 Double_t dlong, Double_t dlat, Double_t dpsi)
1685{
1686 if (!IsPerspective()) return;
1687 if (nsteps<1) return;
1688 Double_t fc = 1./Double_t(nsteps);
1689 Double_t oc[3], od[3], dir[3];
1690 dir[0] = 0;
1691 dir[1] = 0;
1692 dir[2] = 1.;
1693 Int_t i, j;
1694 for (i=0; i<3; i++) {
1695 oc[i] = 0.5*(fRmin[i]+fRmax[i]);
1696 od[i] = 0.5*(fRmax[i]-fRmin[i]);
1697 }
1698 Double_t dox = cov[0]-oc[0];
1699 Double_t doy = cov[1]-oc[1];
1700 Double_t doz = cov[2]-oc[2];
1701
1702 Double_t dd = TMath::Sqrt(dox*dox+doy*doy+doz*doz);
1703 if (dd!=0) {;
1704 dir[0] = dox/dd;
1705 dir[1] = doy/dd;
1706 dir[2] = doz/dd;
1707 }
1708 dd *= fc;
1709 dox = fc*(dx-od[0]);
1710 doy = fc*(dy-od[1]);
1711 doz = fc*(dz-od[2]);
1712 for (i=0; i<nsteps; i++) {
1713 oc[0] += dd*dir[0];
1714 oc[1] += dd*dir[1];
1715 oc[2] += dd*dir[2];
1716 od[0] += dox;
1717 od[1] += doy;
1718 od[2] += doz;
1719 for (j=0; j<3; j++) {
1720 fRmin[j] = oc[j]-od[j];
1721 fRmax[j] = oc[j]+od[j];
1722 }
1724 fLatitude += dlat;
1725 fLongitude += dlong;
1726 fPsi += dpsi;
1728 if (gPad) {
1729 gPad->Modified();
1730 gPad->Update();
1731 }
1732 }
1733}
1734
1735////////////////////////////////////////////////////////////////////////////////
1736/// - 'a' increase scale factor (clip cube borders)
1737/// - 's' decrease scale factor (clip cube borders)
1738
1740{
1741 if (count <= 0) count = 1;
1742 switch (option) {
1743 case '+':
1744 ZoomView();
1745 break;
1746 case '-':
1747 UnzoomView();
1748 break;
1749 case 's':
1750 case 'S':
1751 UnzoomView();
1752 break;
1753 case 'a':
1754 case 'A':
1755 ZoomView();
1756 break;
1757 case 'l':
1758 case 'L':
1759 case 'h':
1760 case 'H':
1761 case 'u':
1762 case 'U':
1763 case 'i':
1764 case 'I':
1765 MoveWindow(option);
1766 break;
1767 case 'j':
1768 case 'J':
1769 ZoomIn();
1770 break;
1771 case 'k':
1772 case 'K':
1773 ZoomOut();
1774 break;
1775 default:
1776 break;
1777 }
1778}
1779
1780////////////////////////////////////////////////////////////////////////////////
1781/// Move view window :
1782/// - l,L - left
1783/// - h,H - right
1784/// - u,U - down
1785/// - i,I - up
1786
1788{
1789 if (!IsPerspective()) return;
1790 Double_t shiftu = 0.1*fUVcoord[2];
1791 Double_t shiftv = 0.1*fUVcoord[3];
1792 switch (option) {
1793 case 'l':
1794 case 'L':
1795 fUVcoord[0] += shiftu;
1796 break;
1797 case 'h':
1798 case 'H':
1799 fUVcoord[0] -= shiftu;
1800 break;
1801 case 'u':
1802 case 'U':
1803 fUVcoord[1] += shiftv;
1804 break;
1805 case 'i':
1806 case 'I':
1807 fUVcoord[1] -= shiftv;
1808 break;
1809 default:
1810 return;
1811 }
1813 if (gPad) {
1814 gPad->Modified();
1815 gPad->Update();
1816 }
1817}
1818
1819////////////////////////////////////////////////////////////////////////////////
1820/// Zoom in.
1821
1823{
1824 if (!IsPerspective()) return;
1825 Double_t extent = GetExtent();
1826 Double_t fc = 0.1;
1827 if (fDview<extent) {
1828 fDview -= fc*extent;
1829 } else {
1830 fDview /= 1.25;
1831 }
1833 if (gPad) {
1834 gPad->Modified();
1835 gPad->Update();
1836 }
1837}
1838
1839////////////////////////////////////////////////////////////////////////////////
1840/// Zoom out.
1841
1843{
1844 if (!IsPerspective()) return;
1845 Double_t extent = GetExtent();
1846 Double_t fc = 0.1;
1847 if (fDview<extent) {
1848 fDview += fc*extent;
1849 } else {
1850 fDview *= 1.25;
1851 }
1853 if (gPad) {
1854 gPad->Modified();
1855 gPad->Update();
1856 }
1857}
1858
1859////////////////////////////////////////////////////////////////////////////////
1860/// Stream an object of class TView3D.
1861
1862void TView3D::Streamer(TBuffer &R__b)
1863{
1864 if (R__b.IsReading()) {
1865 UInt_t R__s, R__c;
1866 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
1867 if (R__v > 1) {
1868 R__b.ReadClassBuffer(TView3D::Class(), this, R__v, R__s, R__c);
1869 return;
1870 }
1871 //====process old versions before automatic schema evolution
1872 //unfortunately we forgot to increment the TView3D version number
1873 //when the class was upgraded to double precision in version 2.25.
1874 //we are forced to use the file version number to recognize old files.
1875 if (R__b.GetParent() && R__b.GetVersionOwner() < 22500) { //old version in single precision
1876 TObject::Streamer(R__b);
1877 TAttLine::Streamer(R__b);
1878 Float_t single, sa[12];
1879 Int_t i;
1880 R__b >> fSystem;
1881 R__b >> single; fLatitude = single;
1882 R__b >> single; fLongitude = single;
1883 R__b >> single; fPsi = single;
1884 R__b.ReadStaticArray(sa); for (i=0;i<12;i++) fTN[i] = sa[i];
1885 R__b.ReadStaticArray(sa); for (i=0;i<12;i++) fTB[i] = sa[i];
1886 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fRmax[i] = sa[i];
1887 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fRmin[i] = sa[i];
1888 R__b.ReadStaticArray(sa); for (i=0;i<12;i++) fTnorm[i] = sa[i];
1889 R__b.ReadStaticArray(sa); for (i=0;i<12;i++) fTback[i] = sa[i];
1890 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fX1[i] = sa[i];
1891 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fX2[i] = sa[i];
1892 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fY1[i] = sa[i];
1893 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fY2[i] = sa[i];
1894 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fZ1[i] = sa[i];
1895 R__b.ReadStaticArray(sa); for (i=0;i<3;i++) fZ2[i] = sa[i];
1896 R__b >> fOutline;
1897 R__b >> fDefaultOutline;
1898 R__b >> fAutoRange;
1899 } else {
1900 TObject::Streamer(R__b);
1901 TAttLine::Streamer(R__b);
1902 R__b >> fLatitude;
1903 R__b >> fLongitude;
1904 R__b >> fPsi;
1905 R__b.ReadStaticArray(fTN);
1906 R__b.ReadStaticArray(fTB);
1907 R__b.ReadStaticArray(fRmax);
1908 R__b.ReadStaticArray(fRmin);
1909 R__b.ReadStaticArray(fTnorm);
1910 R__b.ReadStaticArray(fTback);
1911 R__b.ReadStaticArray(fX1);
1912 R__b.ReadStaticArray(fX2);
1913 R__b.ReadStaticArray(fY1);
1914 R__b.ReadStaticArray(fY2);
1915 R__b.ReadStaticArray(fZ1);
1916 R__b.ReadStaticArray(fZ2);
1917 R__b >> fSystem;
1918 R__b >> fOutline;
1919 R__b >> fDefaultOutline;
1920 R__b >> fAutoRange;
1921 }
1922 //====end of old versions
1923
1924 } else {
1925 R__b.WriteClassBuffer(TView3D::Class(),this);
1926 }
1927}
1928
1929// Shortcuts for menus
1938
1939
@ kMouseMotion
Definition: Buttons.h:23
@ kKeyPress
Definition: Buttons.h:20
@ kButton1Motion
Definition: Buttons.h:20
@ kButton1Up
Definition: Buttons.h:19
@ kButton1Down
Definition: Buttons.h:17
void Class()
Definition: Class.C:29
@ kRotate
Definition: GuiTypes.h:373
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define c(i)
Definition: RSha256.hxx:101
#define s1(x)
Definition: RSha256.hxx:91
#define e(i)
Definition: RSha256.hxx:103
static const double x2[5]
static const double x1[5]
int Int_t
Definition: RtypesCore.h:43
short Version_t
Definition: RtypesCore.h:63
char Char_t
Definition: RtypesCore.h:31
unsigned int UInt_t
Definition: RtypesCore.h:44
const Bool_t kFALSE
Definition: RtypesCore.h:90
double Double_t
Definition: RtypesCore.h:57
float Float_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:89
#define ClassImp(name)
Definition: Rtypes.h:361
float xmin
Definition: THbookFile.cxx:93
float ymin
Definition: THbookFile.cxx:93
float xmax
Definition: THbookFile.cxx:93
float ymax
Definition: THbookFile.cxx:93
#define gROOT
Definition: TROOT.h:406
const Double_t kRad
Definition: TView3D.cxx:35
const Int_t kCARTESIAN
Definition: TView3D.cxx:33
const Int_t kPOLAR
Definition: TView3D.cxx:34
#define gPad
Definition: TVirtualPad.h:287
#define gVirtualX
Definition: TVirtualX.h:338
static struct mg_connection * fc(struct mg_context *ctx)
Definition: civetweb.c:3728
static TAxis3D * ToggleZoom(TVirtualPad *pad=0)
Turn ON / OFF the "Ruler" and "zoom mode" of the TAxis3D object attached to the current pad (if pad =...
Definition: TAxis3D.cxx:765
static TAxis3D * ToggleRulers(TVirtualPad *pad=0)
Turn ON / OFF the "Ruler", TAxis3D object attached to the current pad.
Definition: TAxis3D.cxx:737
Buffer base class used for serializing objects.
Definition: TBuffer.h:42
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
TObject * GetParent() const
Return pointer to parent of this buffer.
Definition: TBuffer.cxx:262
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
virtual Int_t GetVersionOwner() const =0
Bool_t IsReading() const
Definition: TBuffer.h:85
virtual Int_t ReadStaticArray(Bool_t *b)=0
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
virtual void Delete(Option_t *option="")=0
Delete this object.
virtual void Paint(Option_t *option="")
Paint all objects in this collection.
A doubly linked list.
Definition: TList.h:44
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:357
TObject & operator=(const TObject &rhs)
TObject assignment operator.
Definition: TObject.h:283
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:694
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:891
@ kMustCleanup
if object destructor must call RecursiveRemove()
Definition: TObject.h:60
static void DrawOutlineCube(TList *outline, Double_t *rmin, Double_t *rmax)
Draw cube outline with 3d polylines.
The 3D view class.
Definition: TView3D.h:29
virtual void Zoom()
Definition: TView3D.cxx:1936
Double_t fDproj
Definition: TView3D.h:36
Bool_t fDefaultOutline
Definition: TView3D.h:54
virtual void UnzoomView(TVirtualPad *pad=0, Double_t unZoomFactor=1.25)
unZOOM this view.
Definition: TView3D.cxx:1651
Double_t fTB[16]
Definition: TView3D.h:40
TView3D()
Default constructor.
Definition: TView3D.cxx:109
Double_t fUpix
Definition: TView3D.h:37
virtual void Side()
Definition: TView3D.cxx:1933
Double_t fX2[3]
Definition: TView3D.h:47
Bool_t fChanged
Definition: TView3D.h:56
static void AdjustPad(TVirtualPad *pad=0)
Force the current pad to be updated.
Definition: TView3D.cxx:1539
virtual void ZoomView(TVirtualPad *pad=0, Double_t zoomFactor=1.25)
ZOOM this view.
Definition: TView3D.cxx:1660
virtual void SetPerspective()
Set perspective option.
Definition: TView3D.cxx:1315
Double_t fPsi
Definition: TView3D.h:34
virtual ~TView3D()
TView3D default destructor.
Definition: TView3D.cxx:276
virtual void Centered3DImages(TVirtualPad *pad=0)
Move view into the center of the scene.
Definition: TView3D.cxx:1635
virtual Int_t GetSystem()
Definition: TView3D.h:101
virtual void SetParallel()
Set the parallel option (default).
Definition: TView3D.cxx:1304
Double_t fLongitude
Definition: TView3D.h:33
virtual void RotateView(Double_t phi, Double_t theta, TVirtualPad *pad=0)
API to rotate view and adjust the pad provided it the current one.
Definition: TView3D.cxx:1555
Double_t fZ1[3]
Definition: TView3D.h:50
Double_t fUVcoord[4]
Definition: TView3D.h:43
virtual void FindPhiSectors(Int_t iopt, Int_t &kphi, Double_t *aphi, Int_t &iphi1, Int_t &iphi2)
Find critical PHI sectors.
Definition: TView3D.cxx:806
Double_t fZ2[3]
Definition: TView3D.h:51
TView3D & operator=(const TView3D &)
Assignment operator.
Definition: TView3D.cxx:235
Double_t fY1[3]
Definition: TView3D.h:48
virtual void PadRange(Int_t rback)
Set the correct window size for lego and surface plots.
Definition: TView3D.cxx:1158
virtual Int_t GetDistancetoAxis(Int_t axis, Int_t px, Int_t py, Double_t &ratio)
Return distance to axis from point px,py.
Definition: TView3D.cxx:971
virtual void WCtoNDC(const Float_t *pw, Float_t *pn)
Transfer point from world to normalized coordinates.
Definition: TView3D.cxx:1479
Double_t fTback[16]
Definition: TView3D.h:45
Double_t fTnorm[16]
Definition: TView3D.h:44
void ResetView(Double_t longitude, Double_t latitude, Double_t psi, Int_t &irep)
Set view direction (in spherical coordinates).
Definition: TView3D.cxx:1437
virtual void DrawOutlineCube(TList *outline, Double_t *rmin, Double_t *rmax)
Draw the outline of a cube while rotating a 3-d object in the pad.
Definition: TView3D.cxx:620
Double_t fLatitude
Definition: TView3D.h:32
virtual void NormalWCtoNDC(const Float_t *pw, Float_t *pn)
Transfer vector of NORMAL from word to normalized coordinates.
Definition: TView3D.cxx:1100
virtual void MoveFocus(Double_t *center, Double_t dx, Double_t dy, Double_t dz, Int_t nsteps=10, Double_t dlong=0, Double_t dlat=0, Double_t dpsi=0)
Move focus to a different box position and extent in nsteps.
Definition: TView3D.cxx:1683
virtual void DefinePerspectiveView()
Define perspective view.
Definition: TView3D.cxx:400
virtual Double_t GetPsi()
Definition: TView3D.h:92
virtual void AdjustScales(TVirtualPad *pad=0)
Adjust all sides of view in respect of the biggest one.
Definition: TView3D.cxx:1617
virtual void MoveWindow(Char_t option)
Move view window :
Definition: TView3D.cxx:1787
virtual void SetAxisNDC(const Double_t *x1, const Double_t *x2, const Double_t *y1, const Double_t *y2, const Double_t *z1, const Double_t *z2)
Store axis coordinates in the NDC system.
Definition: TView3D.cxx:1237
virtual void SetWindow(Double_t u0, Double_t v0, Double_t du, Double_t dv)
Set viewing window.
Definition: TView3D.cxx:1386
virtual void SideView(TVirtualPad *pad=0)
Set to side view.
Definition: TView3D.cxx:1576
virtual void SetView(Double_t longitude, Double_t latitude, Double_t psi, Int_t &irep)
Set view parameters.
Definition: TView3D.cxx:1397
virtual void FindScope(Double_t *scale, Double_t *center, Int_t &irep)
Find centre of a MIN-MAX scope and scale factors.
Definition: TView3D.cxx:933
virtual void SetRange(const Double_t *min, const Double_t *max)
Set Range function.
Definition: TView3D.cxx:1327
@ kPerspective
Definition: TView3D.h:67
virtual Bool_t IsClippedNDC(Double_t *p) const
Check if point is clipped in perspective view.
Definition: TView3D.cxx:1056
virtual void TopView(TVirtualPad *pad=0)
Set to top view.
Definition: TView3D.cxx:1592
virtual void Top()
Definition: TView3D.cxx:1934
virtual void NDCtoWC(const Float_t *pn, Float_t *pw)
Transfer point from normalized to world coordinates.
Definition: TView3D.cxx:1070
virtual void ResizePad()
Recompute window for perspective view.
Definition: TView3D.cxx:1405
Double_t fY2[3]
Definition: TView3D.h:49
virtual void SetDefaultWindow()
Set default viewing window.
Definition: TView3D.cxx:1252
virtual void Front()
Definition: TView3D.cxx:1931
Double_t fX1[3]
Definition: TView3D.h:46
virtual void GetWindow(Double_t &u0, Double_t &v0, Double_t &du, Double_t &dv) const
Get current window extent.
Definition: TView3D.cxx:1045
virtual void ExecuteRotateView(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TView3D.cxx:642
virtual void FindNormal(Double_t x, Double_t y, Double_t z, Double_t &zn)
Find Z component of NORMAL in normalized coordinates.
Definition: TView3D.cxx:786
Double_t fTN[16]
Definition: TView3D.h:39
virtual void UnZoom()
Definition: TView3D.cxx:1937
virtual Bool_t IsPerspective() const
Definition: TView3D.h:109
Bool_t fAutoRange
Definition: TView3D.h:55
virtual void ZoomMove()
Definition: TView3D.cxx:1935
virtual void ToggleZoom(TVirtualPad *pad=0)
Turn on /off the interactive option to Zoom / Move / Change attributes of 3D axis correspond this vie...
Definition: TView3D.cxx:1609
Double_t fRmin[3]
Definition: TView3D.h:42
Double_t fRmax[3]
Definition: TView3D.h:41
virtual void MoveViewCommand(Char_t chCode, Int_t count=1)
Definition: TView3D.cxx:1739
virtual void ToggleRulers(TVirtualPad *pad=0)
Turn on /off 3D axis.
Definition: TView3D.cxx:1600
virtual void ZoomOut()
Zoom out.
Definition: TView3D.cxx:1842
Double_t fDview
Definition: TView3D.h:35
virtual void Centered()
Definition: TView3D.cxx:1930
virtual void FrontView(TVirtualPad *pad=0)
Set to front view.
Definition: TView3D.cxx:1584
virtual Double_t GetExtent() const
Get maximum view extent.
Definition: TView3D.cxx:1017
TSeqCollection * fOutline
Definition: TView3D.h:53
virtual void DefineViewDirection(const Double_t *s, const Double_t *c, Double_t cosphi, Double_t sinphi, Double_t costhe, Double_t sinthe, Double_t cospsi, Double_t sinpsi, Double_t *tnorm, Double_t *tback)
Define view direction (in spherical coordinates)
Definition: TView3D.cxx:514
virtual void GetRange(Float_t *min, Float_t *max)
Get Range function.
Definition: TView3D.cxx:1029
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TView3D.cxx:628
virtual void ZoomIn()
Zoom in.
Definition: TView3D.cxx:1822
virtual void FindThetaSectors(Int_t iopt, Double_t phi, Int_t &kth, Double_t *ath, Int_t &ith1, Int_t &ith2)
Find critical THETA sectors for given PHI sector.
Definition: TView3D.cxx:872
Int_t fSystem
Definition: TView3D.h:52
virtual void SetOutlineToCube()
This is a function which creates default outline.
Definition: TView3D.cxx:1292
Double_t fVpix
Definition: TView3D.h:38
virtual void AxisVertex(Double_t ang, Double_t *av, Int_t &ix1, Int_t &ix2, Int_t &iy1, Int_t &iy2, Int_t &iz1, Int_t &iz2)
Define axis vertices.
Definition: TView3D.cxx:298
virtual void ShowAxis()
Definition: TView3D.cxx:1932
See TView3D.
Definition: TView.h:25
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition: TVirtualPad.h:51
virtual void Modified(Bool_t flag=1)=0
virtual void AbsCoordinates(Bool_t set)=0
virtual void Update()=0
virtual void SetTheta(Double_t theta=30)=0
virtual void SetPhi(Double_t phi=30)=0
Abstract 3D shapes viewer.
return c1
Definition: legend1.C:41
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
return c2
Definition: legend2.C:14
return c3
Definition: legend3.C:15
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
static constexpr double s
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
Double_t Sqrt(Double_t x)
Definition: TMath.h:681
Double_t Cos(Double_t)
Definition: TMath.h:631
Double_t Sin(Double_t)
Definition: TMath.h:627
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
auto * th2
Definition: textalign.C:17
auto * th1
Definition: textalign.C:13
auto * a
Definition: textangle.C:12