From d0484b81c7ac86b2e65cbeef23e44cf1629c5743 Mon Sep 17 00:00:00 2001 From: Rene Brun <Rene.Brun@cern.ch> Date: Fri, 7 Feb 2003 13:46:48 +0000 Subject: [PATCH] New version from Andrei and Mihaela: - several implementations of method Safety for shapes - few bug fixes in safety computation (phi segmented shapes) - volumes have now an unique ID that you can get with GetNumber() - Mihaela implemented an overlap checker that can be called either for all geometry : TGeoManager::CheckOverlaps(Double_t ovlp=0.1) or just for a given volume. OVLP is the lower overlapping limit. The checker looks for the following errors in the geometry definition : 1 - nodes extruding their mother - can be visualised with TGeoManager::DrawExtrusion(char *mother, char *node) which draws mother in blue and extruding daughter in red 2 - ONLY brothers overlapping to each other - can be visualized with TGeoManager::DrawOverlap(char *volume, char *node1, char *node2) - node1 will appear in red and node2 in blue. These methods can be called only after CheckOverlaps() that prints the list of overlaps. They are all in the context menu of both TGeoVolume and TGeoManager. git-svn-id: http://root.cern.ch/svn/root/trunk@6068 27541ba8-7e3a-0410-8455-c3a389f83636 --- geom/inc/TGeoBBox.h | 3 +- geom/inc/TGeoCone.h | 4 +- geom/inc/TGeoManager.h | 7 +- geom/inc/TGeoPcon.h | 3 +- geom/inc/TGeoShape.h | 4 +- geom/inc/TGeoSphere.h | 3 +- geom/inc/TGeoTube.h | 4 +- geom/inc/TGeoVolume.h | 6 +- geom/inc/TVirtualGeoPainter.h | 8 +- geom/src/TGeoArb8.cxx | 60 ++++--- geom/src/TGeoBBox.cxx | 9 +- geom/src/TGeoCone.cxx | 46 ++++-- geom/src/TGeoManager.cxx | 43 ++++- geom/src/TGeoPara.cxx | 42 +++-- geom/src/TGeoPcon.cxx | 78 ++++++--- geom/src/TGeoPgon.cxx | 192 ++++++++++++++++------ geom/src/TGeoShape.cxx | 20 ++- geom/src/TGeoSphere.cxx | 55 ++++++- geom/src/TGeoTrd1.cxx | 268 ++++++++++++------------------ geom/src/TGeoTrd2.cxx | 280 ++++++++++++-------------------- geom/src/TGeoTube.cxx | 62 ++++--- geom/src/TGeoVolume.cxx | 83 +++++++++- geompainter/inc/TGeoChecker.h | 3 +- geompainter/inc/TGeoPainter.h | 11 ++ geompainter/src/TGeoChecker.cxx | 135 ++++++++++++++- geompainter/src/TGeoPainter.cxx | 112 +++++++++++++ 26 files changed, 1021 insertions(+), 520 deletions(-) diff --git a/geom/inc/TGeoBBox.h b/geom/inc/TGeoBBox.h index 459c3f5b055..d91477b621b 100644 --- a/geom/inc/TGeoBBox.h +++ b/geom/inc/TGeoBBox.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoBBox.h,v 1.4 2002/09/27 16:16:05 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoBBox.h,v 1.5 2003/01/23 14:25:36 brun Exp $ // Author: Andrei Gheata 24/10/01 /************************************************************************* @@ -70,6 +70,7 @@ public: virtual Bool_t IsCylType() const {return kFALSE;} virtual Bool_t IsValidBox() const {return ((fDX<0)||(fDY<0)||(fDZ<0))?kFALSE:kTRUE;} virtual Bool_t IsNullBox() const {return ((fDX==0)&&(fDY==0)&&(fDZ==0))?kTRUE:kFALSE;} + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); diff --git a/geom/inc/TGeoCone.h b/geom/inc/TGeoCone.h index ede0ae5c57c..9346acf4ad8 100644 --- a/geom/inc/TGeoCone.h +++ b/geom/inc/TGeoCone.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoCone.h,v 1.7 2003/01/23 14:25:36 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoCone.h,v 1.8 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -82,6 +82,7 @@ public: virtual void InspectShape() const; virtual Bool_t IsCylType() const {return kTRUE;} + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); @@ -151,6 +152,7 @@ public: Double_t GetPhi1() const {return fPhi1;} Double_t GetPhi2() const {return fPhi2;} virtual void InspectShape() const; + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); diff --git a/geom/inc/TGeoManager.h b/geom/inc/TGeoManager.h index 8a9ac85b567..ffdfed91a13 100644 --- a/geom/inc/TGeoManager.h +++ b/geom/inc/TGeoManager.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoManager.h,v 1.24 2003/01/28 18:12:56 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoManager.h,v 1.25 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 25/10/01 /************************************************************************* @@ -111,7 +111,7 @@ public: Int_t AddMaterial(const TGeoMaterial *material); Int_t AddTransformation(const TGeoMatrix *matrix); Int_t AddShape(const TGeoShape *shape); - Int_t AddVolume(const TGeoVolume *volume); + Int_t AddVolume(TGeoVolume *volume); //--- browsing and tree navigation void Browse(TBrowser *b); virtual Bool_t cd(const char *path=""); // *MENU* @@ -149,6 +149,9 @@ public: //--- geometry checking void CheckGeometry(Option_t *option=""); + void CheckOverlaps(Double_t ovlp=0.1, Option_t *option=""); // *MENU* + void DrawOverlap(const char *mother, const char *node1, const char *node2); // *MENU* + void DrawExtrusion(const char *mother, const char *node); // *MENU* void CheckPoint(Double_t x=0,Double_t y=0, Double_t z=0, Option_t *option=""); // *MENU* void DrawCurrentPoint(Int_t color=2); // *MENU* void DrawPath(const char *path); diff --git a/geom/inc/TGeoPcon.h b/geom/inc/TGeoPcon.h index 33e97a2967b..30d6e9da5fe 100644 --- a/geom/inc/TGeoPcon.h +++ b/geom/inc/TGeoPcon.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoPcon.h,v 1.5 2002/12/03 16:01:39 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoPcon.h,v 1.6 2003/01/23 14:25:36 brun Exp $ // Author: Andrei Gheata 24/10/01 /************************************************************************* @@ -73,6 +73,7 @@ public: virtual TGeoShape *GetMakeRuntimeShape(TGeoShape * /*mother*/) const {return 0;} virtual void InspectShape() const; virtual Bool_t IsCylType() const {return kTRUE;} + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); diff --git a/geom/inc/TGeoShape.h b/geom/inc/TGeoShape.h index 1227357c4e9..a5c5bbc3dd7 100644 --- a/geom/inc/TGeoShape.h +++ b/geom/inc/TGeoShape.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoShape.h,v 1.6 2002/12/06 16:45:03 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoShape.h,v 1.7 2003/01/23 14:25:36 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -113,10 +113,12 @@ public: Bool_t IsValid() const {return !TestBit(kGeoInvalidShape);} virtual Bool_t IsValidBox() const = 0; virtual void InspectShape() const = 0; + virtual void *Make3DBuffer(const TGeoVolume *vol) const = 0; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const = 0; virtual void Paint(Option_t *option) = 0; virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option) = 0; virtual Double_t Safety(Double_t *point, Bool_t in=kTRUE) const = 0; + static Double_t SafetyPhi(Double_t *point, Bool_t in, Double_t c1, Double_t s1, Double_t c2, Double_t s2); virtual void SetDimensions(Double_t *param) = 0; void SetId(Int_t id) {fShapeId = id;} virtual void SetPoints(Double_t *buff) const = 0; diff --git a/geom/inc/TGeoSphere.h b/geom/inc/TGeoSphere.h index 8066e36fc68..2aeb188c46e 100644 --- a/geom/inc/TGeoSphere.h +++ b/geom/inc/TGeoSphere.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoSphere.h,v 1.6 2002/12/06 16:45:03 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoSphere.h,v 1.7 2003/01/23 14:25:36 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -77,6 +77,7 @@ public: virtual void InspectShape() const; virtual Bool_t IsCylType() const {return kFALSE;} Bool_t IsPointInside(Double_t *point, Bool_t checkR=kTRUE, Bool_t checkTh=kTRUE, Bool_t checkPh=kTRUE) const; + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); diff --git a/geom/inc/TGeoTube.h b/geom/inc/TGeoTube.h index 010c14c0713..b91b6ed4055 100644 --- a/geom/inc/TGeoTube.h +++ b/geom/inc/TGeoTube.h @@ -1,4 +1,4 @@ -// @(#)root/base:$Name: $:$Id: TGeoTube.h,v 1.7 2003/01/23 14:25:36 brun Exp $ +// @(#)root/base:$Name: $:$Id: TGeoTube.h,v 1.8 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 24/10/01 /************************************************************************* @@ -67,6 +67,7 @@ public: virtual Double_t GetDz() const {return fDz;} virtual void InspectShape() const; virtual Bool_t IsCylType() const {return kTRUE;} + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); @@ -131,6 +132,7 @@ public: Double_t GetPhi1() const {return fPhi1;} Double_t GetPhi2() const {return fPhi2;} virtual void InspectShape() const; + virtual void *Make3DBuffer(const TGeoVolume *vol) const; virtual void NextCrossing(TGeoParamCurve *c, Double_t *point) const; virtual void Paint(Option_t *option); virtual void PaintNext(TGeoHMatrix *glmat, Option_t *option); diff --git a/geom/inc/TGeoVolume.h b/geom/inc/TGeoVolume.h index 50b6a57e7ec..00b831bf1a8 100644 --- a/geom/inc/TGeoVolume.h +++ b/geom/inc/TGeoVolume.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoVolume.h,v 1.20 2003/01/27 13:16:26 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoVolume.h,v 1.21 2003/01/27 18:04:47 brun Exp $ // Author: Andrei Gheata 30/05/02 /************************************************************************* @@ -94,6 +94,9 @@ public: void ClearShape(); void CleanAll(); void CheckGeometry(Int_t nrays=1, Double_t startx=0, Double_t starty=0, Double_t startz=0) const; + void CheckOverlaps(Double_t ovlp=0.1, Option_t *option="") const; // *MENU* + void DrawExtrusion(const char *node); // *MENU* + void DrawOverlap(const char *node1, const char *node2); // *MENU* Int_t CountNodes(Int_t nlevels=1000) const; // *MENU* Bool_t Contains(Double_t *point) const {return fShape->Contains(point);} Bool_t IsFolder() const; @@ -142,6 +145,7 @@ public: Bool_t IsStyleDefault() const; void InspectMaterial() const; // *MENU* void InspectShape() const {fShape->InspectShape();} // *MENU* + void *Make3DBuffer() const {return fShape->Make3DBuffer(this);} TGeoVolume *MakeCopyVolume(); void MakeCopyNodes(const TGeoVolume *other); Bool_t OptimizeVoxels(); // *MENU* diff --git a/geom/inc/TVirtualGeoPainter.h b/geom/inc/TVirtualGeoPainter.h index a90d1dbd2ce..262b3c8c0e9 100644 --- a/geom/inc/TVirtualGeoPainter.h +++ b/geom/inc/TVirtualGeoPainter.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TVirtualGeoPainter.h,v 1.11 2003/01/06 17:05:43 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TVirtualGeoPainter.h,v 1.12 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 11/01/02 /************************************************************************* @@ -64,6 +64,7 @@ public: virtual void BombTranslation(const Double_t *tr, Double_t *bombtr) = 0; virtual void CheckPoint(Double_t x=0, Double_t y=0, Double_t z=0, Option_t *option="") = 0; virtual void CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, Double_t startz) const = 0; + virtual void CheckOverlaps(const TGeoVolume *vol, Double_t ovlp=0.1, Option_t *option="") const = 0; virtual void DefaultAngles() = 0; virtual void DefaultColors() = 0; virtual Int_t DistanceToPrimitiveVol(TGeoVolume *vol, Int_t px, Int_t py) = 0; @@ -86,6 +87,11 @@ public: virtual TH2F *LegoPlot(Int_t ntheta=60, Double_t themin=0., Double_t themax=180., Int_t nphi=90, Double_t phimin=0., Double_t phimax=360., Double_t rmin=0., Double_t rmax=9999999, Option_t *option="") = 0; + virtual void *MakeBox3DBuffer(const TGeoVolume *vol) = 0; + virtual void *MakeTube3DBuffer(const TGeoVolume *vol) = 0; + virtual void *MakeTubs3DBuffer(const TGeoVolume *vol) = 0; + virtual void *MakePcon3DBuffer(const TGeoVolume *vol) = 0; + virtual void *MakeSphere3DBuffer(const TGeoVolume *vol) = 0; virtual void ModifiedPad() const = 0; virtual void Paint(Option_t *option="") = 0; virtual void PaintBox(TGeoShape *shape, Option_t *option="", TGeoHMatrix *glmat=0) = 0; diff --git a/geom/src/TGeoArb8.cxx b/geom/src/TGeoArb8.cxx index 60e46db076c..57f027eef53 100644 --- a/geom/src/TGeoArb8.cxx +++ b/geom/src/TGeoArb8.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoArb8.cxx,v 1.17 2003/01/23 14:25:36 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoArb8.cxx,v 1.18 2003/01/24 08:38:50 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -768,33 +768,41 @@ Double_t TGeoTrap::Safety(Double_t *point, Bool_t in) const Double_t saf[5]; Double_t norm[3]; // normal to current facette Int_t i; // current facette index - Double_t x0, y0, x1, y1; // lower coordinater - Double_t x0h, y0h, x1h, y1h; // upper coordinates - Double_t dl, dh; // lengths of lower and upper segments - Double_t rl, rh, rlh; - Double_t st, ct, sp, cp; // Sin/Cos of spherical angles of current normal + Double_t x0, y0, z0, x1, y1, z1, x2, y2, z2; + Double_t ax, ay, az, bx, by; + Double_t fn; //---> compute safety for lateral planes for (i=0; i<4; i++) { x0 = fXY[i][0]; y0 = fXY[i][1]; - x1 = fXY[(i+1)%4][0]; - y1 = fXY[(i+1)%4][1]; - x0h = fXY[i+4][0]; - y0h = fXY[i+4][1]; - x1h = fXY[4+((i+1)%4)][0]; - y1h = fXY[4+((i+1)%4)][1]; - dl = TMath::Sqrt((x1-x0)*(x1-x0)+(y1-y0)*(y1-y0)); - dh = TMath::Sqrt((x1h-x0h)*(x1h-x0h)+(y1h-y0h)*(y1h-y0h)); - rl = TMath::Abs(x0*y1-x1*y0)/dl; - rh = TMath::Abs(x0h*y1h-x1h*y0h)/dh; - rlh = TMath::Sqrt((rl-rh)*(rl-rh)+4*fDz*fDz); - st = 2*fDz/rlh; - ct = (rl-rh)/rlh; - sp = (x1-x0)/dl; - cp = (y0-y1)/dl; - norm[0] = st*cp; - norm[1] = st*sp; - norm[2] = ct; + z0 = -fDz; + x1 = fXY[i+4][0]; + y1 = fXY[i+4][1]; + z1 = fDz; + ax = x1-x0; + ay = y1-y0; + az = z1-z0; + x2 = fXY[(i+1)%4][0]; + y2 = fXY[(i+1)%4][1]; + z2 = -fDz; + bx = x2-x0; + by = y2-y0; + if (bx==0 && by==0) { + x2 = fXY[4+((i+1)%4)][0]; + y2 = fXY[4+((i+1)%4)][1]; + z2 = fDz; + bx = x2-x1; + by = y2-y1; + if (bx==0 && by==0) continue; + } + norm[0] = -az*by; + norm[1] = az*bx; + norm[2] = ax*by-ay*bx; + fn = TMath::Sqrt(norm[0]*norm[0]+norm[1]*norm[1]+norm[2]*norm[2]); + if (fn<1E-10) continue; + norm[0] /= fn; + norm[1] /= fn; + norm[2] /= fn; saf[i] = (x0-point[0])*norm[0]+(y0-point[1])*norm[1]+(-fDz-point[2])*norm[2]; if (in) { saf[i]=TMath::Abs(saf[i]); // they should be all positive anyway @@ -802,11 +810,11 @@ Double_t TGeoTrap::Safety(Double_t *point, Bool_t in) const saf[i] = -saf[i]; // only negative values are interesting } } - saf[4] = TMath::Abs(fDz-TMath::Abs(point[2])); + saf[4] = fDz-TMath::Abs(point[2]); if (in) { safe = saf[TMath::LocMin(5, saf)]; } else { - if (TMath::Abs(point[2]) <= fDz) saf[4]=-saf[4]; + saf[4]=-saf[4]; safe = saf[TMath::LocMax(5, saf)]; } return safe; diff --git a/geom/src/TGeoBBox.cxx b/geom/src/TGeoBBox.cxx index 583796768a2..fa24aab0f11 100644 --- a/geom/src/TGeoBBox.cxx +++ b/geom/src/TGeoBBox.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoBBox.cxx,v 1.15 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoBBox.cxx,v 1.16 2003/01/27 13:16:26 brun Exp $ // Author: Andrei Gheata 24/10/01 // Contains() and DistToIn/Out() implemented by Mihaela Gheata @@ -348,6 +348,13 @@ void TGeoBBox::InspectShape() const printf(" origin: x=%11.5f y=%11.5f z=%11.5f\n", fOrigin[0], fOrigin[1], fOrigin[2]); } //----------------------------------------------------------------------------- +void *TGeoBBox::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeBox3DBuffer(vol); +} +//----------------------------------------------------------------------------- void TGeoBBox::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const { // computes next intersection point of curve c with this shape diff --git a/geom/src/TGeoCone.cxx b/geom/src/TGeoCone.cxx index 7dfc039414c..b627379ac3c 100644 --- a/geom/src/TGeoCone.cxx +++ b/geom/src/TGeoCone.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoCone.cxx,v 1.14 2003/01/27 13:16:26 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoCone.cxx,v 1.15 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoCone::Contains() and DistToOut() implemented by Mihaela Gheata @@ -487,6 +487,13 @@ void TGeoCone::InspectShape() const TGeoBBox::InspectShape(); } //----------------------------------------------------------------------------- +void *TGeoCone::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeTube3DBuffer(vol); +} +//----------------------------------------------------------------------------- void TGeoCone::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const { // computes next intersection point of curve c with this shape @@ -526,7 +533,7 @@ Double_t TGeoCone::Safety(Double_t *point, Bool_t in) const Double_t rin = tg1*point[2]+ro1; Double_t rout = tg2*point[2]+ro2; saf[0] = fDz-TMath::Abs(point[2]); - saf[1] = (ro1>0)?((r-rin)*cr1):(-kBig); + saf[1] = (ro1>0)?((r-rin)*cr1):kBig; saf[2] = (rout-r)*cr2; if (in) return saf[TMath::LocMin(3,saf)]; for (Int_t i=0; i<3; i++) saf[i]=-saf[i]; @@ -561,7 +568,7 @@ Double_t TGeoCone::SafetyS(Double_t *point, Bool_t in, Double_t dz, Double_t rmi default: saf[0] = dz-TMath::Abs(point[2]); } - saf[1] = (ro1>0)?((r-rin)*cr1):(-kBig); + saf[1] = (ro1>0)?((r-rin)*cr1):kBig; saf[2] = (rout-r)*cr2; if (in) return saf[TMath::LocMin(3,saf)]; for (Int_t i=0; i<3; i++) saf[i]=-saf[i]; @@ -1262,6 +1269,13 @@ void TGeoConeSeg::InspectShape() const TGeoBBox::InspectShape(); } //----------------------------------------------------------------------------- +void *TGeoConeSeg::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeTubs3DBuffer(vol); +} +//----------------------------------------------------------------------------- void TGeoConeSeg::Paint(Option_t *option) { // paint this shape according to option @@ -1291,7 +1305,7 @@ Double_t TGeoConeSeg::Safety(Double_t *point, Bool_t in) const // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - Double_t saf[5]; + Double_t saf[4]; Double_t phi1 = fPhi1*kDegRad; Double_t phi2 = fPhi2*kDegRad; Double_t c1 = TMath::Cos(phi1); @@ -1310,21 +1324,20 @@ Double_t TGeoConeSeg::Safety(Double_t *point, Bool_t in) const Double_t rout = tg2*point[2]+ro2; saf[0] = fDz-TMath::Abs(point[2]); // positive if inside - saf[1] = (ro1>1E-10)?((r-rin)*cr1):kBig; + saf[1] = (r-rin)*cr1; saf[2] = (rout-r)*cr2; - saf[3] = -point[0]*s1 + point[1]*c1; - saf[4] = point[0]*s2 - point[1]*c2; + saf[3] = TGeoShape::SafetyPhi(point, in, c1,s1,c2,s2); - if (in) return saf[TMath::LocMin(5,saf)]; - for (Int_t i=0; i<5; i++) saf[i]=-saf[i]; - return saf[TMath::LocMax(5,saf)]; + if (in) return saf[TMath::LocMin(4,saf)]; + for (Int_t i=0; i<4; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(4,saf)]; } //----------------------------------------------------------------------------- Double_t TGeoConeSeg::SafetyS(Double_t *point, Bool_t in, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Int_t skipz) { // Static method to compute the closest distance from given point to this shape. - Double_t saf[5]; + Double_t saf[4]; Double_t ro1 = 0.5*(rmin1+rmin2); Double_t tg1 = 0.5*(rmin2-rmin1)/dz; Double_t cr1 = 1./TMath::Sqrt(1.+tg1*tg1); @@ -1348,14 +1361,13 @@ Double_t TGeoConeSeg::SafetyS(Double_t *point, Bool_t in, Double_t dz, Double_t default: saf[0] = dz-TMath::Abs(point[2]); } - saf[1] = (ro1>1E-10)?((r-rin)*cr1):kBig; + saf[1] = (r-rin)*cr1; saf[2] = (rout-r)*cr2; - saf[3] = -point[0]*s1 + point[1]*c1; - saf[4] = point[0]*s2 - point[1]*c2; + saf[3] = TGeoShape::SafetyPhi(point,in,c1,s1,c2,s2); - if (in) return saf[TMath::LocMin(5,saf)]; - for (Int_t i=0; i<5; i++) saf[i]=-saf[i]; - return saf[TMath::LocMax(5,saf)]; + if (in) return saf[TMath::LocMin(4,saf)]; + for (Int_t i=0; i<4; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(4,saf)]; } //----------------------------------------------------------------------------- diff --git a/geom/src/TGeoManager.cxx b/geom/src/TGeoManager.cxx index 47e5641adf7..7aaa9e43d3c 100644 --- a/geom/src/TGeoManager.cxx +++ b/geom/src/TGeoManager.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoManager.cxx,v 1.40 2003/01/29 10:42:20 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoManager.cxx,v 1.41 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 25/10/01 /************************************************************************* @@ -603,7 +603,7 @@ Int_t TGeoManager::AddShape(const TGeoShape *shape) return index; } //_____________________________________________________________________________ -Int_t TGeoManager::AddVolume(const TGeoVolume *volume) +Int_t TGeoManager::AddVolume(TGeoVolume *volume) { // Add a volume to the list. Returns index of the volume in list. if (!volume) { @@ -614,6 +614,7 @@ Int_t TGeoManager::AddVolume(const TGeoVolume *volume) if (volume->IsRunTime()) list = fGVolumes; Int_t index = list->GetSize(); list->Add((TGeoVolume*)volume); + volume->SetNumber(index); return index; } //_____________________________________________________________________________ @@ -2675,6 +2676,7 @@ void TGeoManager::SetNsegments(Int_t nseg) { // Set number of segments for approximating circles in drawing. GetGeomPainter(); + if (fNsegments==nseg) return; if (nseg>2) fNsegments = nseg; if (fPainter) fPainter->SetNsegments(nseg); } @@ -2824,6 +2826,43 @@ void TGeoManager::CheckGeometry(Option_t * /*option*/) fTopNode->CheckShapes(); } +//_____________________________________________________________________________ +void TGeoManager::CheckOverlaps(Double_t ovlp, Option_t * option) +{ +// Check all geometry for illegal overlaps within a limit OVLP. + TIter next(fVolumes); + TGeoVolume *vol; + while ((vol=(TGeoVolume*)next())) { + if (!vol->GetNdaughters() || vol->GetFinder()) continue; + vol->CheckOverlaps(ovlp, option); + } +} + +//_____________________________________________________________________________ +void TGeoManager::DrawExtrusion(const char *mother, const char *node) +{ +// Draw togeather only a given volume and one daughter node, as given by CheckOverlaps(). +// This method offer a visual validation of a declared extrusion (daughter not fully +// contained by its mother) + TGeoVolume *vol = (TGeoVolume*)fVolumes->FindObject(mother); + if (!vol) { + Error("DrawExtrusion", "volume %s not found", mother); + return; + } + vol->DrawExtrusion(node); +} + +//_____________________________________________________________________________ +void TGeoManager::DrawOverlap(const char *mother, const char *node1, const char *node2) +{ +// Draw togeather only 2 possible overlapping daughters of a given volume. + TGeoVolume *vol = (TGeoVolume*)fVolumes->FindObject(mother); + if (!vol) { + Error("DrawExtrusion", "volume %s not found", mother); + return; + } + vol->DrawOverlap(node1, node2); +} //_____________________________________________________________________________ void TGeoManager::UpdateCurrentPosition(Double_t * /*nextpoint*/) { diff --git a/geom/src/TGeoPara.cxx b/geom/src/TGeoPara.cxx index 7e27a0195cc..e1f74dfb793 100644 --- a/geom/src/TGeoPara.cxx +++ b/geom/src/TGeoPara.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoPara.cxx,v 1.11 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoPara.cxx,v 1.12 2003/01/27 13:16:26 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoPara::Contains() implemented by Mihaela Gheata @@ -143,6 +143,12 @@ Bool_t TGeoPara::Contains(Double_t *point) const Double_t TGeoPara::DistToOut(Double_t *point, Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const { // compute distance from inside point to surface of the para + if (iact<3 && safe) { + // compute safety + *safe = Safety(point, kTRUE); + if (iact==0) return kBig; + if (iact==1 && step<*safe) return kBig; + } Double_t saf[6]; Double_t snxt = kBig; // distance from point to higher Z face @@ -156,25 +162,12 @@ Double_t TGeoPara::DistToOut(Double_t *point, Double_t *dir, Int_t iact, Double_ saf[2] = fY-yt; // distance from point to lower Y face saf[3] = -fY-yt; - // cos of angle YZ - Double_t cty = 1.0/TMath::Sqrt(1.0+fTyz*fTyz); - // distance from point to center axis on X Double_t xt = point[0]-fTxz*point[2]-fTxy*yt; // distance from point to higher X face saf[0] = fX-xt; // distance from point to lower X face saf[1] = -fX-xt; - // cos of angle XZ - Double_t ctx = 1.0/TMath::Sqrt(1.0+fTxy*fTxy+fTxz*fTxz); - if (iact<3 && safe) { - // compute safety - *safe = TMath::Min(saf[0]*ctx, -saf[1]*ctx); - *safe = TMath::Min(*safe, TMath::Min(saf[2]*cty, -saf[3]*cty)); - *safe = TMath::Min(*safe, TMath::Min(saf[4], -saf[5])); - if (iact==0) return kBig; - if (iact==1 && step<*safe) return kBig; - } Double_t sn1, sn2, sn3; sn1 = sn2 = sn3 = kBig; if (dir[2]!=0) sn3=saf[4]/dir[2]; @@ -403,11 +396,28 @@ void TGeoPara::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const // computes next intersection point of curve c with this shape } //----------------------------------------------------------------------------- -Double_t TGeoPara::Safety(Double_t * /*point*/, Bool_t /*in*/) const +Double_t TGeoPara::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - return kBig; + Double_t saf[3]; + // distance from point to higher Z face + saf[0] = fZ-TMath::Abs(point[2]); // Z + + Double_t yt = point[1]-fTyz*point[2]; + saf[1] = fY-TMath::Abs(yt); // Y + // cos of angle YZ + Double_t cty = 1.0/TMath::Sqrt(1.0+fTyz*fTyz); + + Double_t xt = point[0]-fTxz*point[2]-fTxy*yt; + saf[2] = fX-TMath::Abs(xt); // X + // cos of angle XZ + Double_t ctx = 1.0/TMath::Sqrt(1.0+fTxy*fTxy+fTxz*fTxz); + saf[2] *= ctx; + saf[1] *= cty; + if (in) return saf[TMath::LocMin(3,saf)]; + for (Int_t i=0; i<3; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(3,saf)]; } //----------------------------------------------------------------------------- void TGeoPara::SetDimensions(Double_t *param) diff --git a/geom/src/TGeoPcon.cxx b/geom/src/TGeoPcon.cxx index 92b86a08490..544213286e7 100644 --- a/geom/src/TGeoPcon.cxx +++ b/geom/src/TGeoPcon.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoPcon.cxx,v 1.14 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoPcon.cxx,v 1.15 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 24/10/01 // TGeoPcon::Contains() implemented by Mihaela Gheata @@ -53,6 +53,7 @@ TGeoPcon::TGeoPcon(Double_t phi, Double_t dphi, Int_t nz) // Default constructor SetBit(TGeoShape::kGeoPcon); fPhi1 = phi; + if (fPhi1<0) fPhi1+=360.; fDphi = dphi; fNz = nz; fRmin = new Double_t [nz]; @@ -66,6 +67,7 @@ TGeoPcon::TGeoPcon(const char *name, Double_t phi, Double_t dphi, Int_t nz) // Default constructor SetBit(TGeoShape::kGeoPcon); fPhi1 = phi; + if (fPhi1<0) fPhi1+=360.; fDphi = dphi; fNz = nz; fRmin = new Double_t [nz]; @@ -109,7 +111,6 @@ void TGeoPcon::ComputeBBox() rmax = fRmax[TMath::LocMax(fNz, fRmax)]; Double_t phi1 = fPhi1; Double_t phi2 = phi1 + fDphi; - if (phi2 > 360) phi2-=360; Double_t xc[4]; Double_t yc[4]; @@ -129,19 +130,15 @@ void TGeoPcon::ComputeBBox() Double_t ddp = -phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) xmax = rmax; ddp = 90-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) ymax = rmax; ddp = 180-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) xmin = -rmax; ddp = 270-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) ymin = -rmax; fOrigin[0] = (xmax+xmin)/2; fOrigin[1] = (ymax+ymin)/2; @@ -165,7 +162,7 @@ Bool_t TGeoPcon::Contains(Double_t *point) const while ((izh-izl)>1) { if (point[2] > fZ[izt]) izl = izt; else izh = izt; - izt = (izl+izh)/2; + izt = (izl+izh)>>1; } // the point is in the section bounded by izl and izh Z planes @@ -182,12 +179,14 @@ Bool_t TGeoPcon::Contains(Double_t *point) const } if ((r2<rmin*rmin) || (r2>rmax*rmax)) return kFALSE; // now check phi - Double_t phi = fPhi1+0.5*fDphi; - if ((point[1]!=0.0) || (point[0] != 0.0)) - phi = TMath::ATan2(point[1], point[0]) * kRadDeg; - if (phi < fPhi1) phi+=360.0; - if ((phi<fPhi1) || ((phi-fPhi1)>fDphi)) return kFALSE; - return kTRUE; + if (fDphi==360) return kTRUE; + if (r2<1E-10) return kTRUE; + Double_t phi = TMath::ATan2(point[1], point[0]) * kRadDeg; + if (phi < 0) phi+=360.0; + Double_t ddp = phi-fPhi1; + if (ddp<0) ddp+=360.; + if (ddp<=fDphi) return kTRUE; + return kFALSE; } //----------------------------------------------------------------------------- Int_t TGeoPcon::DistancetoPrimitive(Int_t px, Int_t py) @@ -520,6 +519,14 @@ void TGeoPcon::InspectShape() const printf(" plane %i: z=%11.5f Rmin=%11.5f Rmax=%11.5f\n", ipl, fZ[ipl], fRmin[ipl], fRmax[ipl]); TGeoBBox::InspectShape(); } +//----------------------------------------------------------------------------- +void *TGeoPcon::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakePcon3DBuffer(vol); +} + //----------------------------------------------------------------------------- void TGeoPcon::Paint(Option_t *option) { @@ -556,6 +563,8 @@ Double_t TGeoPcon::Safety(Double_t *point, Bool_t in) const Bool_t is_tube, is_seg; Double_t phi1=0, phi2=0, c1=0, s1=0, c2=0, s2=0; Int_t skipz; + Double_t saf[2]; + saf[0] = saf[1] = kBig; if (in) { //---> point is inside pcon Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]); @@ -563,23 +572,38 @@ Double_t TGeoPcon::Safety(Double_t *point, Bool_t in) const if (ipl<0) return 0; // point on first Z boundary dz = 0.5*(fZ[ipl+1]-fZ[ipl]); if (dz<1E-10) return 0; - Double_t dzb = TMath::Min(point[2]-fZ[ipl], fZ[ipl+1]-point[2]); - if (dzb<1E-5) return dzb; // point close to Z boundary skipz = 3; // skip z checks if (ipl==0) { - if (fNz>2 && fZ[1]==fZ[2]) skipz=0; // do not skip z plane check - else skipz=2; // skip upper z check - } else if (ipl==fNz-2) { - if (fNz>2 && fZ[ipl]==fZ[ipl-1]) skipz=0; - else skipz=1; // skip lower z check - } else { + saf[0] = point[2]-fZ[0]; + if (saf[0]<1E-4) return saf[0]; + } + if (ipl==fNz-2) { + saf[1] = fZ[fNz-1]-point[2]; + if (saf[1]<1E-4) return saf[1]; + } + if (ipl>1) { if (fZ[ipl]==fZ[ipl-1]) { - skipz=2; - if (fZ[ipl+1]==fZ[ipl+2]) skipz=0; - } else if (fZ[ipl+1]==fZ[ipl+2]) { - skipz=1; + if (fRmin[ipl]<fRmin[ipl-1] || fRmax[ipl]>fRmax[ipl-1]) { + saf[0] = point[2]-fZ[ipl]; + if (saf[0]<1E-4) return saf[0]; + } } - } + } + if (ipl<fNz-3) { + if (fZ[ipl+1]==fZ[ipl+2]) { + if (fRmin[ipl+1]<fRmin[ipl+2] || fRmax[ipl+1]>fRmax[ipl+2]) { + saf[1] = fZ[ipl+1]-point[2]; + if (saf[1]<1E-4) return saf[1]; + } + } + } + if (saf[0]<1E10) { + if (saf[1]<1E10) skipz=0; // check both Z planes + else skipz=2; // skip upper Z + } else { + if (saf[1]<1E10) skipz=1; // skip lower Z + else skipz=3; // skip both Z planes + } //---> Check shape type memcpy(ptnew, point, 3*sizeof(Double_t)); ptnew[2] -= 0.5*(fZ[ipl]+fZ[ipl+1]); @@ -673,6 +697,7 @@ Double_t TGeoPcon::Safety(Double_t *point, Bool_t in) const rmax2 = fRmax[ipnew+1]; is_tube = ((rmin1==rmin2) && (rmax1==rmax2))?kTRUE:kFALSE; ptnew[2] = point[2] - 0.5*(fZ[ipnew]+fZ[ipnew+1]); + dz = 0.5*(fZ[ipl+1]-fZ[ipl]); if (is_seg) { if (is_tube) safdown = TGeoTubeSeg::SafetyS(ptnew,in,rmin1,rmax1, dz,c1,s1,c2,s2,skipz); else safdown = TGeoConeSeg::SafetyS(ptnew,in,dz,rmin1,rmax1,rmin2,rmax2,c1,s1,c2,s2,skipz); @@ -698,6 +723,7 @@ Double_t TGeoPcon::Safety(Double_t *point, Bool_t in) const rmax2 = fRmax[ipnew+1]; is_tube = ((rmin1==rmin2) && (rmax1==rmax2))?kTRUE:kFALSE; ptnew[2] = point[2] - 0.5*(fZ[ipnew]+fZ[ipnew+1]); + dz = 0.5*(fZ[ipl+1]-fZ[ipl]); if (is_seg) { if (is_tube) safup = TGeoTubeSeg::SafetyS(ptnew,in,rmin1,rmax1, dz,c1,s1,c2,s2,skipz); else safup = TGeoConeSeg::SafetyS(ptnew,in,dz,rmin1,rmax1,rmin2,rmax2,c1,s1,c2,s2,skipz); diff --git a/geom/src/TGeoPgon.cxx b/geom/src/TGeoPgon.cxx index 17664c29434..a1293724646 100644 --- a/geom/src/TGeoPgon.cxx +++ b/geom/src/TGeoPgon.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoPgon.cxx,v 1.19 2003/01/27 13:16:26 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoPgon.cxx,v 1.20 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoPgon::Contains() implemented by Mihaela Gheata @@ -97,7 +97,6 @@ void TGeoPgon::ComputeBBox() rmax = rmax/TMath::Cos(0.5*divphi*kDegRad); Double_t phi1 = fPhi1; Double_t phi2 = phi1 + fDphi; - if (phi2 > 360) phi2-=360; Double_t xc[4]; Double_t yc[4]; @@ -117,19 +116,15 @@ void TGeoPgon::ComputeBBox() Double_t ddp = -phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) xmax = rmax; ddp = 90-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) ymax = rmax; ddp = 180-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) xmin = -rmax; ddp = 270-phi1; if (ddp<0) ddp+= 360; - if (ddp>360) ddp-=360; if (ddp<=fDphi) ymin = -rmax; fOrigin[0] = (xmax+xmin)/2; fOrigin[1] = (ymax+ymin)/2; @@ -143,46 +138,33 @@ Bool_t TGeoPgon::Contains(Double_t *point) const { // test if point is inside this shape // check total z range - if ((point[2]<fZ[0]) || (point[2]>fZ[fNz-1])) return kFALSE; - // find smallest Rmin and largest Rmax - Double_t rmin = fRmin[0]; - Double_t rmax = fRmax[0]; - for (Int_t i=1; i<fNz; i++) { - if (fRmin[i] < rmin) rmin = fRmin[i]; - if (fRmax[i] > rmax) rmax = fRmax[i]; - } - // check R against rmin - Double_t r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]); - if (r < rmin) return kFALSE; + if (point[2]<fZ[0]) return kFALSE; + if (point[2]>fZ[fNz-1]) return kFALSE; Double_t divphi = fDphi/fNedges; - // find the radius of the outscribed circle - rmax = rmax/TMath::Cos(0.5*divphi*kDegRad); - // check R against rmax - if (r > rmax) return kFALSE; // now check phi Double_t phi = TMath::ATan2(point[1], point[0])*kRadDeg; - if (phi < fPhi1) phi += 360.0; - if ((phi<fPhi1) || ((phi-fPhi1)>fDphi)) return kFALSE; + if (phi < 0) phi += 360.0; + Double_t ddp = phi-fPhi1; + if (ddp<0) ddp+=360.; + if (ddp>fDphi) return kFALSE; // now find phi division - Int_t ipsec = (Int_t)TMath::Min((phi-fPhi1)/divphi+1., (Double_t)fNedges); - Double_t ph0 = (fPhi1+divphi*(ipsec-0.5))*kDegRad; + Int_t ipsec = TMath::Min(Int_t(ddp/divphi), fNedges-1); + Double_t ph0 = (fPhi1+divphi*(ipsec+0.5))*kDegRad; // now check projected distance - r = point[0]*TMath::Cos(ph0) + point[1]*TMath::Sin(ph0); + Double_t r = point[0]*TMath::Cos(ph0) + point[1]*TMath::Sin(ph0); // find in which Z section the point is in - Int_t izl = 0; - Int_t izh = fNz-1; - Int_t izt = 0; - while ((izh-izl)>1) { - izt = (izl+izh)/2; - if (point[2] < fZ[izt]) izh = izt; - else izl=izt; - } + Int_t iz = TMath::BinarySearch(fNz, fZ, point[2]); + if (iz==fNz-1) { + if (r<fRmin[iz]) return kFALSE; + if (r>fRmax[iz]) return kFALSE; + return kTRUE; + } // now compute rmin and rmax and test the value of r - Double_t dzrat = (point[2]-fZ[izl])/(fZ[izh]-fZ[izl]); - rmin = fRmin[izl]+dzrat*(fRmin[izh]-fRmin[izl]); + Double_t dzrat = (point[2]-fZ[iz])/(fZ[iz+1]-fZ[iz]); + Double_t rmin = fRmin[iz]+dzrat*(fRmin[iz+1]-fRmin[iz]); // is the point inside the 'hole' at the center of the volume ? if (r < rmin) return kFALSE; - rmax = fRmax[izl]+dzrat*(fRmax[izh]-fRmax[izl]); + Double_t rmax = fRmax[iz]+dzrat*(fRmax[iz+1]-fRmax[iz]); if (r > rmax) return kFALSE; return kTRUE; @@ -556,15 +538,24 @@ Double_t TGeoPgon::DistToOutSect(Double_t *point, Double_t *dir, Int_t &iz, Int_ return dmin; } //----------------------------------------------------------------------------- -Double_t TGeoPgon::DistToOut(Double_t *point, Double_t *dir, Int_t /*iact*/, Double_t /*step*/, Double_t * /*safe*/) const +Double_t TGeoPgon::DistToOut(Double_t *point, Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const { // compute distance from inside point to surface of the polygone // first find out in which Z section the point is in + if (iact<3 && safe) { + *safe = Safety(point, kTRUE); + if (iact==0) return kBig; + if (iact==1 && step<*safe) return kBig; + } Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]); - if ((ipl==(fNz-1)) || (ipl<0)) { + if (ipl==fNz-1) { + ipl--; + if (dir[2]>=0) return 1E-6; + } + if (ipl<0) { // point out - Warning("DistToOut", "point is outside Z range"); - return kBig; + ipl++; + if (dir[2]<=0) return 1E-6; } // Double_t dz = 0.5*(fZ[ipl+1]-fZ[ipl]); @@ -585,6 +576,11 @@ Double_t TGeoPgon::DistToIn(Double_t *point, Double_t *dir, Int_t iact, Double_t { // compute distance from outside point to surface of the polygone // first find in which segment we are + if (iact<3 && safe) { + *safe = Safety(point, kFALSE); + if (iact==0) return kBig; + if (iact==1 && step<*safe) return kBig; + } Double_t pt[3]; Double_t eps = 0; memcpy(&pt[0], point, 3*sizeof(Double_t)); @@ -686,11 +682,6 @@ Double_t TGeoPgon::DistToIn(Double_t *point, Double_t *dir, Int_t iact, Double_t Double_t dph2=(bits & kInphi)?(phi2-phi):(phi-phi2); saf[4]=r*TMath::Sin(dph1*kDegRad); saf[5]=r*TMath::Sin(dph2*kDegRad); - if ((iact<3) && safe) { - *safe = saf[TMath::LocMax(6, &saf[0])]; -// if ((iact==1) && (*safe>step)) return kBig; - if (iact==0) return kBig; - } // compute distance to boundary if (!cross) return kBig; Double_t snxt=DistToInSect(&pt[0], dir, ifirst, ipsec, bits, &saf[0]); @@ -983,11 +974,116 @@ void TGeoPgon::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const // computes next intersection point of curve c with this shape } //----------------------------------------------------------------------------- -Double_t TGeoPgon::Safety(Double_t * /*point*/, Bool_t /*in*/) const +Double_t TGeoPgon::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - return kBig; + Double_t saf[5]; + Double_t safe, dz; + Int_t i; + Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]); + for (i=0; i<5; i++) saf[i]=kBig; + if (in) { + //---> first locate Z segment and compute Z safety + if (ipl==(fNz-1)) return 0; + if (ipl<0) return 0; + dz = fZ[ipl+1]-fZ[ipl]; + if (dz<1E-10) return 0; + if (ipl==0) { + saf[0] = point[2]-fZ[0]; + if (saf[0]<1E-4) return saf[0]; + } + if (ipl==fNz-2) { + saf[1] = fZ[fNz-1]-point[2]; + if (saf[1]<1E-4) return saf[1]; + } + if (ipl>1) { + if (fZ[ipl]==fZ[ipl-1]) { + if (fRmin[ipl]<fRmin[ipl-1] || fRmax[ipl]>fRmax[ipl-1]) { + saf[0] = point[2]-fZ[ipl]; + if (saf[0]<1E-4) return saf[0]; + } + } + } + if (ipl<fNz-3) { + if (fZ[ipl+1]==fZ[ipl+2]) { + if (fRmin[ipl+1]<fRmin[ipl+2] || fRmax[ipl+1]>fRmax[ipl+2]) { + saf[1] = fZ[ipl+1]-point[2]; + if (saf[1]<1E-4) return saf[1]; + } + } + } + } else { + if (ipl>=0 && ipl<fNz-1) { + dz = fZ[ipl+1]-fZ[ipl]; + if (dz==0) { + ipl++; + dz = fZ[ipl+1]-fZ[ipl]; + } + if (ipl>1) { + if (fZ[ipl]==fZ[ipl-1]) { + if (fRmin[ipl]>fRmin[ipl-1] || fRmax[ipl]<fRmax[ipl-1]) { + saf[0] = point[2]-fZ[ipl]; + if (saf[0]<1E-4) return saf[0]; + } + } + } + if (ipl<fNz-3) { + if (fZ[ipl+1]==fZ[ipl+2]) { + if (fRmin[ipl+1]>fRmin[ipl+2] || fRmax[ipl+1]<fRmax[ipl+2]) { + saf[1] = fZ[ipl+1]-point[2]; + if (saf[1]<1E-4) return saf[1]; + } + } + } + } else { + if (ipl<0) { + ipl=0; + saf[0] = fZ[0]-point[2]; + } else { + ipl=fNz-2; + saf[1] = point[2]-fZ[fNz-1]; + } + dz = fZ[ipl+1]-fZ[ipl]; + } + } + //---> compute phi safety + if (fDphi<360) { + Double_t phi1 = fPhi1*kDegRad; + Double_t phi2 = (fPhi1+fDphi)*kDegRad; + Double_t c1 = TMath::Cos(phi1); + Double_t s1 = TMath::Sin(phi1); + Double_t c2 = TMath::Cos(phi2); + Double_t s2 = TMath::Sin(phi2); + saf[2] = TGeoShape::SafetyPhi(point,in,c1,s1,c2,s2); + } + + //---> locate phi and compute R safety + Double_t divphi = fDphi/fNedges; + Double_t phi = TMath::ATan2(point[1], point[0])*kRadDeg; + if (phi<0) phi+=360.; + Double_t ddp = phi-fPhi1; + if (ddp<0) ddp+=360.; + Int_t ipsec = Int_t(ddp/divphi); + Double_t ph0 = (fPhi1+divphi*(ipsec+0.5))*kDegRad; + // compute projected distance + Double_t r, rsum, rpgon, ta, calf; + r = point[0]*TMath::Cos(ph0)+point[1]*TMath::Sin(ph0); + rsum = fRmin[ipl]+fRmin[ipl+1]; + if (rsum>1E-10) { + ta = (fRmin[ipl+1]-fRmin[ipl])/dz; + calf = 1./TMath::Sqrt(1+ta*ta); + rpgon = fRmin[ipl] + (point[2]-fZ[ipl])*ta; + saf[3] = (r-rpgon)*calf; + } + ta = (fRmax[ipl+1]-fRmax[ipl])/dz; + calf = 1./TMath::Sqrt(1+ta*ta); + rpgon = fRmax[ipl] + (point[2]-fZ[ipl])*ta; + saf[4] = (rpgon-r)*calf; + if (in) return saf[TMath::LocMin(5,saf)]; + for (i=0; i<5; i++) saf[i]=-saf[i]; + safe = saf[TMath::LocMax(5,saf)]; + return safe; } //----------------------------------------------------------------------------- void TGeoPgon::SetDimensions(Double_t *param) diff --git a/geom/src/TGeoShape.cxx b/geom/src/TGeoShape.cxx index 33ecbb965ef..3dd6fdc13f9 100644 --- a/geom/src/TGeoShape.cxx +++ b/geom/src/TGeoShape.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoShape.cxx,v 1.6 2002/10/08 16:17:49 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoShape.cxx,v 1.7 2003/01/06 17:05:44 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -345,3 +345,21 @@ Int_t TGeoShape::GetVertexNumber(Bool_t vx, Bool_t vy, Bool_t vz) if (!vx) return imax; return imin; } +//----------------------------------------------------------------------------- +Double_t TGeoShape::SafetyPhi(Double_t *point, Bool_t in, Double_t c1, Double_t s1, Double_t c2, Double_t s2) +{ +// Static method to compute safety w.r.t a phi corner defined by cosines/sines +// of the angles phi1, phi2. + Double_t saf1 = kBig; + Double_t saf2 = kBig; + if (point[0]*c1+point[1]*s1 >= 0) saf1 = -point[0]*s1 + point[1]*c1; + if (point[0]*c2+point[1]*s2 >= 0) saf2 = point[0]*s2 - point[1]*c2; + if (in) { + if (saf1<0) saf1=kBig; + if (saf2<0) saf2=kBig; + return TMath::Min(saf1,saf2); + } + if (saf1<0 && saf2<0) return TMath::Max(saf1,saf2); + return TMath::Min(saf1,saf2); +} + diff --git a/geom/src/TGeoSphere.cxx b/geom/src/TGeoSphere.cxx index 8e250b411f1..057a3e0f32f 100644 --- a/geom/src/TGeoSphere.cxx +++ b/geom/src/TGeoSphere.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoSphere.cxx,v 1.13 2003/01/27 13:16:26 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoSphere.cxx,v 1.14 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoSphere::Contains() DistToIn/Out() implemented by Mihaela Gheata @@ -737,17 +737,66 @@ void TGeoSphere::PaintNext(TGeoHMatrix *glmat, Option_t *option) if (!painter) return; painter->PaintSphere(this, option, glmat); } + +//----------------------------------------------------------------------------- +void *TGeoSphere::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeSphere3DBuffer(vol); +} + //----------------------------------------------------------------------------- void TGeoSphere::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const { // computes next intersection point of curve c with this shape } + //----------------------------------------------------------------------------- -Double_t TGeoSphere::Safety(Double_t * /*point*/, Bool_t /*in*/) const +Double_t TGeoSphere::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - return kBig; + Double_t rxy2 = point[0]*point[0]+point[1]*point[1]; + Double_t rxy = TMath::Sqrt(rxy2); + Double_t r2 = rxy2+point[2]*point[2]; + Double_t r=TMath::Sqrt(r2); + Bool_t rzero=kFALSE; + if (r<=1E-20) rzero=kTRUE; + //localize theta + Double_t phi=0;; + Double_t th=0.; + if (TestBit(kGeoThetaSeg) && (!rzero)) { + th = TMath::ACos(point[2]/r)*kRadDeg; + } + //localize phi + if (TestBit(kGeoPhiSeg)) { + phi=TMath::ATan2(point[1], point[0])*kRadDeg; + if (phi<0) phi+=360.; + } + Double_t saf[6]; + saf[0]=(fRmin==0 && !TestBit(kGeoThetaSeg) && !TestBit(kGeoPhiSeg))?kBig:r-fRmin; + saf[1]=fRmax-r; + saf[2]=saf[3]=saf[4]=saf[5]= kBig; + if (TestBit(kGeoThetaSeg)) { + if (fTheta1>0) { + saf[2] = r*TMath::Sin((th-fTheta1)*kDegRad); + } + if (fTheta2<180) { + saf[3] = r*TMath::Sin((fTheta2-th)*kDegRad); + } + } + if (TestBit(kGeoPhiSeg)) { + Double_t dph1=phi-fPhi1; + if (dph1<0) dph1+=360.; + if (dph1<=90.) saf[4]=rxy*TMath::Sin(dph1*kDegRad); + Double_t dph2=fPhi2-phi; + if (dph2<0) dph2+=360.; + if (dph2<=90.) saf[5]=rxy*TMath::Sin(dph2*kDegRad); + } + if (in) return saf[TMath::LocMin(6,saf)]; + for (Int_t i=0; i<6; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(6, saf)]; } //----------------------------------------------------------------------------- void TGeoSphere::SetSphDimensions(Double_t rmin, Double_t rmax, Double_t theta1, diff --git a/geom/src/TGeoTrd1.cxx b/geom/src/TGeoTrd1.cxx index c24af529ecb..af8443da019 100644 --- a/geom/src/TGeoTrd1.cxx +++ b/geom/src/TGeoTrd1.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoTrd1.cxx,v 1.13 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoTrd1.cxx,v 1.14 2003/01/27 13:16:26 brun Exp $ // Author: Andrei Gheata 24/10/01 // TGeoTrd1::Contains() and DistToOut() implemented by Mihaela Gheata @@ -120,70 +120,42 @@ Double_t TGeoTrd1::DistToOut(Double_t *point, Double_t *dir, Int_t iact, Double_ // compute distance from inside point to surface of the trd1 Double_t snxt = kBig; - - Double_t saf[6]; - //--- Compute safety first - // check Z facettes - saf[0] = point[2]+fDz; - saf[1] = fDz-point[2]; - Double_t fx = 0.5*(fDx1-fDx2)/fDz; - Double_t calf = 1./TMath::Sqrt(1.0+fx*fx); - Double_t salf = calf*fx; - Double_t s,cn; - // check X facettes - Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; - saf[2] = (distx+point[0])*calf; - saf[3] = (distx-point[0])*calf; - // check Y facettes - saf[4] = point[1]+fDy; - saf[5] = fDy-point[1]; if (iact<3 && safe) { // compute safe distance - *safe = saf[TMath::LocMin(6, saf)]; + *safe = Safety(point, kTRUE); if (iact==0) return kBig; if (iact==1 && step<*safe) return kBig; } + + //--- Compute safety first + Double_t fx = 0.5*(fDx1-fDx2)/fDz; + Double_t cn; + Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; //--- Compute distance to this shape // first check if Z facettes are crossed - cn = -dir[2]; - gGeoManager->SetNormalChecked(TMath::Abs(cn)); - if (cn>0) { - snxt = saf[0]/cn; - } else { - snxt = -saf[1]/cn; - } + Double_t dist[3]; + for (Int_t i=0; i<3; i++) dist[i]=kBig; + if (dir[2]<0) { + dist[0]=-(point[2]+fDz)/dir[2]; + } else if (dir[2]>0) { + dist[0]=(fDz-point[2])/dir[2]; + } // now check X facettes - cn = -calf*dir[0]+salf*dir[2]; - if (cn>0) { - s = saf[2]/cn; - if (s<snxt) { - snxt = s; - gGeoManager->SetNormalChecked(cn); - } - } - cn = calf*dir[0]+salf*dir[2]; + cn = -dir[0]+fx*dir[2]; + if (cn>0) dist[1] = (point[0]+distx)/cn; + + cn = dir[0]+fx*dir[2]; if (cn>0) { - s = saf[3]/cn; - if (s<snxt) { - snxt = s; - gGeoManager->SetNormalChecked(cn); - } + Double_t s = (distx-point[0])/cn; + if (s<dist[1]) dist[1] = s; } // now check Y facettes - cn = -dir[1]; - if (cn>0) { - s = saf[4]/cn; - if (s<snxt) { - gGeoManager->SetNormalChecked(cn); - return s; - } - } else { - s = -saf[5]/cn; - if (s<snxt) { - gGeoManager->SetNormalChecked(-cn); - return s; - } - } + if (dir[1]<0) { + dist[2]=-(point[1]+fDy)/dir[1]; + } else if (dir[1]>0) { + dist[2]=(fDy-point[1])/dir[1]; + } + snxt = dist[TMath::LocMin(3,dist)]; return snxt; } //----------------------------------------------------------------------------- @@ -252,149 +224,94 @@ Double_t TGeoTrd1::DistToIn(Double_t *point, Double_t *dir, Int_t iact, Double_t { // compute distance from outside point to surface of the trd1 Double_t snxt = kBig; - // find a visible face - Double_t ptnew[3]; - Double_t saf[6]; - memset(saf, 0, 6*sizeof(Double_t)); - //--- Compute safety first - // check visibility of Z facettes - if (point[2]<-fDz) { - saf[0] = -fDz-point[2]; - } else { - if (point[2]>fDz) { - saf[1] = point[2]-fDz; - } - } - Double_t fx = 0.5*(fDx1-fDx2)/fDz; - Double_t calf = 1./TMath::Sqrt(1.0+fx*fx); - Double_t salf = calf*fx; - Double_t cn; - // check visibility of X faces - Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; - if (point[0]<-distx) { - saf[2] = (-point[0]-distx)*calf; - } - if (point[0]>distx) { - saf[3] = (point[0]-distx)*calf; - } - // check visibility of Y facettes - if (point[1]<-fDy) { - saf[4] = -fDy-point[1]; - } else { - if (point[1]>fDy) { - saf[5] = point[1]-fDy; - } - } if (iact<3 && safe) { // compute safe distance - *safe = saf[TMath::LocMax(6, saf)]; + *safe = Safety(point, kFALSE); if (iact==0) return kBig; if (iact==1 && step<*safe) return kBig; } + // find a visible face + Double_t xnew,ynew,znew; + Double_t fx = 0.5*(fDx1-fDx2)/fDz; + Double_t cn; + Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; + //--- Compute distance to this shape // first check if Z facettes are crossed - if (saf[0]>0) { + if (point[2]<-fDz) { cn = -dir[2]; if (cn<0) { - snxt = -saf[0]/cn; + snxt = -(fDz+point[2])/cn; // find extrapolated X and Y - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < fDx1) { - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy) { - // bottom Z facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < fDx1) { + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy) return snxt; } } - } else { - if (saf[1]>0) { - cn = dir[2]; - if (cn<0) { - snxt = -saf[1]/cn; - // find extrapolated X and Y - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < fDx2) { - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy) { - // top Z facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } - } - } - } - } + } else if (point[2]>fDz) { + cn = dir[2]; + if (cn<0) { + snxt = -(fDz-point[2])/cn; + // find extrapolated X and Y + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < fDx2) { + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy) return snxt; + } + } + } // check if X facettes are crossed - if (saf[2]>0) { - cn = -calf*dir[0]+salf*dir[2]; + if (point[0]<-distx) { + cn = -dir[0]+fx*dir[2]; if (cn<0) { - snxt = -saf[2]/cn; + snxt = (point[0]+distx)/cn; // find extrapolated Y and Z - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy) { - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - // lower X facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy) { + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) return snxt; } } } - if (saf[3]>0) { - cn = calf*dir[0]+salf*dir[2]; + if (point[0]>distx) { + cn = dir[0]+fx*dir[2]; if (cn<0) { - snxt = -saf[3]/cn; + snxt = (distx-point[0])/cn; // find extrapolated Y and Z - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy) { - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - // lower X facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy) { + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) return snxt; } } } // finally check Y facettes - if (saf[4]>0) { + if (point[1]<-fDy) { cn = -dir[1]; if (cn<0) { - snxt = -saf[4]/cn; + snxt = (point[1]+fDy)/cn; // find extrapolated X and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - ptnew[0] = point[0]+snxt*dir[0]; - distx = 0.5*(fDx1+fDx2)-fx*ptnew[2]; - if (TMath::Abs(ptnew[0]) < distx) { - // lower Y facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + xnew = point[0]+snxt*dir[0]; + Double_t dx = 0.5*(fDx1+fDx2)-fx*znew; + if (TMath::Abs(xnew) < dx) return snxt; } } - } else { - if (saf[5]>0) { - cn = dir[1]; - if (cn<0) { - snxt = -saf[5]/cn; - // find extrapolated X and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - ptnew[0] = point[0]+snxt*dir[0]; - distx = 0.5*(fDx1+fDx2)-fx*ptnew[2]; - if (TMath::Abs(ptnew[0]) < distx) { - // higher Y facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } - } + } else if (point[1]>0) { + cn = dir[1]; + if (cn<0) { + snxt = (fDy-point[1])/cn; + // find extrapolated X and Z + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + xnew = point[0]+snxt*dir[0]; + Double_t dx = 0.5*(fDx1+fDx2)-fx*znew; + if (TMath::Abs(xnew) < dx) return snxt; } } - } + } return kBig; } //----------------------------------------------------------------------------- @@ -531,12 +448,27 @@ void TGeoTrd1::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const // computes next intersection point of curve c with this shape } //----------------------------------------------------------------------------- -Double_t TGeoTrd1::Safety(Double_t * /*point*/, Bool_t /*in*/) const +Double_t TGeoTrd1::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - return kBig; + Double_t saf[3]; + //--- Compute safety first + // check Z facettes + saf[0] = fDz-TMath::Abs(point[2]); + Double_t fx = 0.5*(fDx1-fDx2)/fDz; + Double_t calf = 1./TMath::Sqrt(1.0+fx*fx); + // check X facettes + Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; + if (distx<0) saf[1]=kBig; + else saf[1]=(distx-TMath::Abs(point[0]))*calf; + // check Y facettes + saf[2] = fDy-TMath::Abs(point[1]); + if (in) return saf[TMath::LocMin(3,saf)]; + for (Int_t i=0; i<3; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(3,saf)]; } + //----------------------------------------------------------------------------- void TGeoTrd1::SetDimensions(Double_t *param) { diff --git a/geom/src/TGeoTrd2.cxx b/geom/src/TGeoTrd2.cxx index fd2e55a21db..7587202de98 100644 --- a/geom/src/TGeoTrd2.cxx +++ b/geom/src/TGeoTrd2.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoTrd2.cxx,v 1.12 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoTrd2.cxx,v 1.13 2003/01/27 13:16:26 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoTrd2::Contains() and DistToOut() implemented by Mihaela Gheata @@ -122,77 +122,46 @@ Double_t TGeoTrd2::DistToOut(Double_t *point, Double_t *dir, Int_t iact, Double_ { // compute distance from inside point to surface of the trd2 Double_t snxt = kBig; - - Double_t saf[6]; - //--- Compute safety first - // check Z facettes - saf[0] = point[2]+fDz; - saf[1] = fDz-point[2]; - Double_t fx = 0.5*(fDx1-fDx2)/fDz; - Double_t fy = 0.5*(fDy1-fDy2)/fDz; - Double_t calfx = 1./TMath::Sqrt(1.0+fx*fx); - Double_t salfx = calfx*fx; - Double_t calfy = 1./TMath::Sqrt(1.0+fy*fy); - Double_t salfy = calfy*fy; - Double_t s,cn; - // check X facettes - Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; - Double_t disty = 0.5*(fDy1+fDy2)-fy*point[2]; - saf[2] = (distx+point[0])*calfx; - saf[3] = (distx-point[0])*calfx; - // check Y facettes - saf[4] = (disty+point[1])*calfy; - saf[5] = (disty-point[1])*calfy; if (iact<3 && safe) { // compute safe distance - *safe = saf[TMath::LocMin(6, saf)]; + *safe = Safety(point, kTRUE); if (iact==0) return kBig; if (iact==1 && step<*safe) return kBig; } + + Double_t fx = 0.5*(fDx1-fDx2)/fDz; + Double_t fy = 0.5*(fDy1-fDy2)/fDz; + Double_t cn; + + Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; + Double_t disty = 0.5*(fDy1+fDy2)-fy*point[2]; + //--- Compute distance to this shape // first check if Z facettes are crossed - cn = -dir[2]; - if (cn>0) { - gGeoManager->SetNormalChecked(cn); - snxt = saf[0]/cn; - } else { - gGeoManager->SetNormalChecked(-cn); - snxt = -saf[1]/cn; - } + Double_t dist[3]; + for (Int_t i=0; i<3; i++) dist[i]=kBig; + if (dir[2]<0) { + dist[0]=-(point[2]+fDz)/dir[2]; + } else if (dir[2]>0) { + dist[0]=(fDz-point[2])/dir[2]; + } // now check X facettes - cn = -calfx*dir[0]+salfx*dir[2]; - if (cn>0) { - s = saf[2]/cn; - if (s<snxt) { - snxt = s; - gGeoManager->SetNormalChecked(cn); - } - } - cn = calfx*dir[0]+salfx*dir[2]; + cn = -dir[0]+fx*dir[2]; + if (cn>0) dist[1] = (point[0]+distx)/cn; + cn = dir[0]+fx*dir[2]; if (cn>0) { - s = saf[3]/cn; - if (s<snxt) { - snxt = s; - gGeoManager->SetNormalChecked(cn); - } + Double_t s = (distx-point[0])/cn; + if (s<dist[1]) dist[1] = s; } // now check Y facettes - cn = -calfy*dir[1]+salfy*dir[2]; + cn = -dir[1]+fy*dir[2]; + if (cn>0) dist[2] = (point[1]+disty)/cn; + cn = dir[1]+fy*dir[2]; if (cn>0) { - s = saf[4]/cn; - if (s<snxt) { - snxt = s; - gGeoManager->SetNormalChecked(cn); - } - } - cn = calfy*dir[1]+salfy*dir[2]; - if (cn>0) { - s = saf[5]/cn; - if (s<snxt) { - gGeoManager->SetNormalChecked(cn); - return s; - } + Double_t s = (disty-point[1])/cn; + if (s<dist[2]) dist[2] = s; } + snxt = dist[TMath::LocMin(3,dist)]; return snxt; } //----------------------------------------------------------------------------- @@ -200,151 +169,96 @@ Double_t TGeoTrd2::DistToIn(Double_t *point, Double_t *dir, Int_t iact, Double_t { // compute distance from outside point to surface of the trd2 Double_t snxt = kBig; + if (iact<3 && safe) { + // compute safe distance + *safe = Safety(point, kFALSE); + if (iact==0) return kBig; + if (iact==1 && step<*safe) return kBig; + } // find a visible face - Double_t ptnew[3]; - Double_t saf[6]; - memset(saf, 0, 6*sizeof(Double_t)); - //--- Compute safety first - // check visibility of Z facettes - if (point[2]<-fDz) { - saf[0] = -fDz-point[2]; - } else { - if (point[2]>fDz) { - saf[1] = point[2]-fDz; - } - } + Double_t xnew,ynew,znew; Double_t fx = 0.5*(fDx1-fDx2)/fDz; - Double_t calfx = 1./TMath::Sqrt(1.0+fx*fx); - Double_t salfx = calfx*fx; Double_t fy = 0.5*(fDy1-fDy2)/fDz; - Double_t calfy = 1./TMath::Sqrt(1.0+fy*fy); - Double_t salfy = calfy*fy; Double_t cn; // check visibility of X faces Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; Double_t disty = 0.5*(fDy1+fDy2)-fy*point[2]; - if (point[0]<-distx) { - saf[2] = (-point[0]-distx)*calfx; - } - if (point[0]>distx) { - saf[3] = (point[0]-distx)*calfx; - } - // check visibility of Y facettes - if (point[1]<-disty) { - saf[4] = (-point[1]-disty)*calfy; - } - if (point[1]>disty) { - saf[5] = (point[1]-disty)*calfy; - } - - if (iact<3 && safe) { - // compute safe distance - *safe = saf[TMath::LocMax(6, saf)]; - if (iact==0) return kBig; - if (iact==1 && step<*safe) return kBig; - } //--- Compute distance to this shape // first check if Z facettes are crossed - if (saf[0]>0) { + if (point[2]<-fDz) { cn = -dir[2]; if (cn<0) { - snxt = -saf[0]/cn; + snxt = -(fDz+point[2])/cn; // find extrapolated X and Y - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < fDx1) { - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy1) { - // bottom Z facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < fDx1) { + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy1) return snxt; } } - } else { - if (saf[1]>0) { - cn = dir[2]; - if (cn<0) { - snxt = -saf[1]/cn; - // find extrapolated X and Y - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < fDx2) { - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < fDy2) { - // top Z facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } - } - } - } - } + } else if (point[2]>fDz) { + cn = dir[2]; + if (cn<0) { + snxt = -(fDz-point[2])/cn; + // find extrapolated X and Y + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < fDx2) { + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < fDy2) return snxt; + } + } + } // check if X facettes are crossed - if (saf[2]>0) { - cn = -calfx*dir[0]+salfx*dir[2]; + if (point[0]<-distx) { + cn = -dir[0]+fx*dir[2]; if (cn<0) { - snxt = -saf[2]/cn; + snxt = (point[0]+distx)/cn; // find extrapolated Y and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - disty = 0.5*(fDy1+fDy2)-fy*ptnew[2]; - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < disty) { - // lower X facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + Double_t dy = 0.5*(fDy1+fDy2)-fy*znew; + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < dy) return snxt; } } } - if (saf[3]>0) { - cn = calfx*dir[0]+salfx*dir[2]; + if (point[0]>distx) { + cn = dir[0]+fx*dir[2]; if (cn<0) { - snxt = -saf[3]/cn; + snxt = (distx-point[0])/cn; // find extrapolated Y and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - disty = 0.5*(fDy1+fDy2)-fy*ptnew[2]; - ptnew[1] = point[1]+snxt*dir[1]; - if (TMath::Abs(ptnew[1]) < disty) { - // upper X facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + Double_t dy = 0.5*(fDy1+fDy2)-fy*znew; + ynew = point[1]+snxt*dir[1]; + if (TMath::Abs(ynew) < dy) return snxt; } } } // finally check Y facettes - if (saf[4]>0) { - cn = -calfy*dir[1]+salfy*dir[2]; + if (point[1]<-disty) { + cn = -dir[1]+fy*dir[2]; if (cn<0) { - snxt = -saf[4]/cn; + snxt = (point[1]+disty)/cn; // find extrapolated X and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - distx = 0.5*(fDx1+fDx2)-fx*ptnew[2]; - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < distx) { - // lower Y facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + Double_t dx = 0.5*(fDx1+fDx2)-fx*znew; + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < dx) return snxt; } } } - if (saf[5]>0) { - cn = calfy*dir[1]+salfy*dir[2]; + if (point[1]>disty) { + cn = dir[1]+fy*dir[2]; if (cn<0) { - snxt = -saf[5]/cn; + snxt = (disty-point[1])/cn; // find extrapolated X and Z - ptnew[2] = point[2]+snxt*dir[2]; - if (TMath::Abs(ptnew[2]) < fDz) { - distx = 0.5*(fDx1+fDx2)-fx*ptnew[2]; - ptnew[0] = point[0]+snxt*dir[0]; - if (TMath::Abs(ptnew[0]) < distx) { - // upper Y facette is crossed - gGeoManager->SetNormalChecked(-cn); - return snxt; - } + znew = point[2]+snxt*dir[2]; + if (TMath::Abs(znew) < fDz) { + Double_t dx = 0.5*(fDx1+fDx2)-fx*znew; + xnew = point[0]+snxt*dir[0]; + if (TMath::Abs(xnew) < dx) return snxt; } } } @@ -540,11 +454,31 @@ void TGeoTrd2::NextCrossing(TGeoParamCurve * /*c*/, Double_t * /*point*/) const // computes next intersection point of curve c with this shape } //----------------------------------------------------------------------------- -Double_t TGeoTrd2::Safety(Double_t * /*point*/, Bool_t /*in*/) const +Double_t TGeoTrd2::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - return kBig; + Double_t saf[3]; + //--- Compute safety first + // check Z facettes + saf[0] = fDz-TMath::Abs(point[2]); + Double_t fx = 0.5*(fDx1-fDx2)/fDz; + Double_t calf = 1./TMath::Sqrt(1.0+fx*fx); + // check X facettes + Double_t distx = 0.5*(fDx1+fDx2)-fx*point[2]; + if (distx<0) saf[1]=kBig; + else saf[1]=(distx-TMath::Abs(point[0]))*calf; + + Double_t fy = 0.5*(fDy1-fDy2)/fDz; + calf = 1./TMath::Sqrt(1.0+fy*fy); + // check Y facettes + distx = 0.5*(fDy1+fDy2)-fy*point[2]; + if (distx<0) saf[2]=kBig; + else saf[2]=(distx-TMath::Abs(point[1]))*calf; + + if (in) return saf[TMath::LocMin(3,saf)]; + for (Int_t i=0; i<3; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(3,saf)]; } //----------------------------------------------------------------------------- void TGeoTrd2::SetDimensions(Double_t *param) diff --git a/geom/src/TGeoTube.cxx b/geom/src/TGeoTube.cxx index 2995aa4e448..b3eda192a4d 100644 --- a/geom/src/TGeoTube.cxx +++ b/geom/src/TGeoTube.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoTube.cxx,v 1.14 2003/01/24 08:38:50 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoTube.cxx,v 1.15 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 24/10/01 // TGeoTube::Contains() and DistToOut/In() implemented by Mihaela Gheata @@ -458,6 +458,14 @@ void TGeoTube::InspectShape() const printf(" dz = %11.5f\n", fDz); TGeoBBox::InspectShape(); } +//----------------------------------------------------------------------------- +void *TGeoTube::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeTube3DBuffer(vol); +} + //----------------------------------------------------------------------------- void TGeoTube::Paint(Option_t *option) { @@ -1127,6 +1135,13 @@ void TGeoTubeSeg::InspectShape() const TGeoBBox::InspectShape(); } //----------------------------------------------------------------------------- +void *TGeoTubeSeg::Make3DBuffer(const TGeoVolume *vol) const +{ + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) return 0; + return painter->MakeTubs3DBuffer(vol); +} +//----------------------------------------------------------------------------- void TGeoTubeSeg::Paint(Option_t *option) { // paint this shape according to option @@ -1154,7 +1169,7 @@ Double_t TGeoTubeSeg::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - Double_t saf[5]; + Double_t saf[4]; Double_t rsq = point[0]*point[0]+point[1]*point[1]; Double_t r = TMath::Sqrt(rsq); Double_t phi1 = fPhi1*kDegRad; @@ -1165,21 +1180,22 @@ Double_t TGeoTubeSeg::Safety(Double_t *point, Bool_t in) const Double_t s2 = TMath::Sin(phi2); saf[0] = fDz-TMath::Abs(point[2]); // positive if inside - saf[1] = (fRmin>1E-10)?(r-fRmin):kBig; + saf[1] = r-fRmin; saf[2] = fRmax-r; - saf[3] = -point[0]*s1 + point[1]*c1; - saf[4] = point[0]*s2 - point[1]*c2; + saf[3] = TGeoShape::SafetyPhi(point,in,c1,s1,c2,s2); + + if (in) return saf[TMath::LocMin(4,saf)]; - if (in) return saf[TMath::LocMin(5,saf)]; - for (Int_t i=0; i<5; i++) saf[i]=-saf[i]; - return saf[TMath::LocMax(5,saf)]; + for (Int_t i=0; i<4; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(4,saf)]; } + //----------------------------------------------------------------------------- Double_t TGeoTubeSeg::SafetyS(Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Int_t skipz) { // Static method to compute the closest distance from given point to this shape. - Double_t saf[5]; + Double_t saf[4]; Double_t rsq = point[0]*point[0]+point[1]*point[1]; Double_t r = TMath::Sqrt(rsq); @@ -1196,15 +1212,16 @@ Double_t TGeoTubeSeg::SafetyS(Double_t *point, Bool_t in, Double_t rmin, Double_ default: saf[0] = dz-TMath::Abs(point[2]); } - saf[1] = (rmin>1E-10)?(r-rmin):kBig; + saf[1] = r-rmin; saf[2] = rmax-r; - saf[3] = -point[0]*s1 + point[1]*c1; - saf[4] = point[0]*s2 - point[1]*c2; + saf[3] = TGeoShape::SafetyPhi(point,in,c1,s1,c2,s2); - if (in) return saf[TMath::LocMin(5,saf)]; - for (Int_t i=0; i<5; i++) saf[i]=-saf[i]; - return saf[TMath::LocMax(5,saf)]; + if (in) return saf[TMath::LocMin(4,saf)]; + + for (Int_t i=0; i<4; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(4,saf)]; } + //----------------------------------------------------------------------------- void TGeoTubeSeg::SetTubsDimensions(Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2) @@ -1819,7 +1836,7 @@ Double_t TGeoCtub::Safety(Double_t *point, Bool_t in) const { // computes the closest distance from given point to this shape, according // to option. The matching point on the shape is stored in spoint. - Double_t saf[6]; + Double_t saf[5]; Double_t rsq = point[0]*point[0]+point[1]*point[1]; Double_t r = TMath::Sqrt(rsq); Bool_t isseg = kTRUE; @@ -1827,7 +1844,7 @@ Double_t TGeoCtub::Safety(Double_t *point, Bool_t in) const saf[0] = -point[0]*fNlow[0] - point[1]*fNlow[1] - (fDz+point[2])*fNlow[2]; saf[1] = -point[0]*fNhigh[0] - point[1]*fNhigh[1] + (fDz-point[2])*fNhigh[2]; - saf[2] = (fRmin>1E-10)?(r-fRmin):kBig; + saf[2] = (fRmin<1E-10 && !isseg)?kBig:(r-fRmin); saf[3] = fRmax-r; if (isseg) { Double_t phi1 = fPhi1*kDegRad; @@ -1837,15 +1854,14 @@ Double_t TGeoCtub::Safety(Double_t *point, Bool_t in) const Double_t c2 = TMath::Cos(phi2); Double_t s2 = TMath::Sin(phi2); - saf[4] = -point[0]*s1 + point[1]*c1; - saf[5] = point[0]*s2 - point[1]*c2; + saf[4] = TGeoShape::SafetyPhi(point, in, c1,s1,c2,s2); } else { - saf[4] = saf[5] = kBig; + saf[4] = kBig; } - if (in) return saf[TMath::LocMin(6,saf)]; - for (Int_t i=0; i<6; i++) saf[i]=-saf[i]; - return saf[TMath::LocMax(6,saf)]; + if (in) return saf[TMath::LocMin(5,saf)]; + for (Int_t i=0; i<5; i++) saf[i]=-saf[i]; + return saf[TMath::LocMax(5,saf)]; } //----------------------------------------------------------------------------- void TGeoCtub::SetCtubDimensions(Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2, diff --git a/geom/src/TGeoVolume.cxx b/geom/src/TGeoVolume.cxx index 70b3f660a78..d35b979a47b 100644 --- a/geom/src/TGeoVolume.cxx +++ b/geom/src/TGeoVolume.cxx @@ -1,6 +1,6 @@ -// @(#)root/geom:$Name: $:$Id: TGeoVolume.cxx,v 1.23 2003/01/27 18:04:47 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoVolume.cxx,v 1.24 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 30/05/02 -// Divide() implemented by Mihaela Gheata +// Divide(), CheckOverlaps() implemented by Mihaela Gheata /************************************************************************* * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. * @@ -179,6 +179,19 @@ void TGeoVolume::CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, Do // if (old_vol) gGeoManager->SetTopVolume(old_vol); } +//_____________________________________________________________________________ +void TGeoVolume::CheckOverlaps(Double_t ovlp, Option_t *option) const +{ +// Overlap checking tool. Check for illegal overlaps within a limit OVLP. + TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); + if (!painter) { + Error("CheckOverlaps", "Could not instanciate painter"); + return; + } + gGeoManager->SetNsegments(80); + painter->CheckOverlaps(this, ovlp, option); +} + //_____________________________________________________________________________ void TGeoVolume::CleanAll() { @@ -499,6 +512,72 @@ void TGeoVolume::DrawOnly(Option_t *option) painter->DrawOnly(option); } +//_____________________________________________________________________________ +void TGeoVolume::DrawExtrusion(const char *node) +{ +// Draw togeather only a given volume and one daughter node, as given by CheckOverlaps(). +// This method offer a visual validation of a declared extrusion (daughter not fully +// contained by its mother) + TGeoNode *daughter = GetNode(node); + if (!daughter) { + Error("DrawExtrusion", "node %s not found", node); + return; + } + Int_t nd = GetNdaughters(); + TGeoNode *current; + for (Int_t i=0; i<nd; i++) { + current = GetNode(i); + if (current==daughter) { + current->SetVisibility(kTRUE); + current->GetVolume()->SetVisibility(kTRUE); + current->GetVolume()->SetLineColor(2); + } else { + current->SetVisibility(kFALSE); + } + } + SetVisibility(kTRUE); + SetLineColor(4); + gGeoManager->SetTopVisible(); + gGeoManager->SetVisLevel(1); + gGeoManager->SetVisOption(0); + Draw(); +} + +//_____________________________________________________________________________ +void TGeoVolume::DrawOverlap(const char *node1, const char *node2) +{ +// Draw togeather only 2 possible overlapping daughters of a given volume. + TGeoNode *daughter1 = GetNode(node1); + if (!daughter1) { + Error("DrawOverlap", "node %s not found", node1); + return; + } + TGeoNode *daughter2 = GetNode(node2); + if (!daughter2) { + Error("DrawOverlap", "node %s not found", node2); + return; + } + Int_t nd = GetNdaughters(); + gGeoManager->SetTopVisible(kFALSE); + gGeoManager->SetVisLevel(1); + gGeoManager->SetVisOption(0); + TGeoNode *current; + for (Int_t i=0; i<nd; i++) { + current = GetNode(i); + if (current==daughter1 || current==daughter2) { + current->SetVisibility(kTRUE); + current->GetVolume()->SetVisibility(kTRUE); + if (current==daughter1) + current->GetVolume()->SetLineColor(2); + else + current->GetVolume()->SetLineColor(4); + } else { + current->SetVisibility(kFALSE); + } + } + Draw(); +} + //_____________________________________________________________________________ Bool_t TGeoVolume::OptimizeVoxels() { diff --git a/geompainter/inc/TGeoChecker.h b/geompainter/inc/TGeoChecker.h index c88695eba98..94ca26cbf0b 100644 --- a/geompainter/inc/TGeoChecker.h +++ b/geompainter/inc/TGeoChecker.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoChecker.h,v 1.5 2002/11/20 08:55:10 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoChecker.h,v 1.6 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 01/11/01 /************************************************************************* @@ -46,6 +46,7 @@ public: virtual ~TGeoChecker(); // methods void CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, Double_t startz) const; + void CheckOverlaps(const TGeoVolume *vol, Double_t ovlp=0.1, Option_t *option="") const; void CheckPoint(Double_t x=0, Double_t y=0, Double_t z=0, Option_t *option=""); Double_t CheckVoxels(TGeoVolume *vol, TGeoVoxelFinder *voxels, Double_t *xyz, Int_t npoints); void CreateTree(const char *treename, const char *filename); diff --git a/geompainter/inc/TGeoPainter.h b/geompainter/inc/TGeoPainter.h index e04626ee864..0ecf26447ec 100644 --- a/geompainter/inc/TGeoPainter.h +++ b/geompainter/inc/TGeoPainter.h @@ -29,6 +29,11 @@ #include "TGeoManager.h" #endif +typedef struct _x3d_points_ { + Int_t numPoints; + Double_t *points; // x0, y0, z0, x1, y1, z1, ... +} X3DPoints; + class TGeoHMatrix; class TGeoChecker; class TH2F; @@ -59,6 +64,7 @@ public: virtual void BombTranslation(const Double_t *tr, Double_t *bombtr); virtual void CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, Double_t startz) const; virtual void CheckPoint(Double_t x=0, Double_t y=0, Double_t z=0, Option_t *option=""); + virtual void CheckOverlaps(const TGeoVolume *vol, Double_t ovlp=0.1, Option_t *option="") const; virtual void DefaultAngles(); virtual void DefaultColors(); virtual Int_t DistanceToPrimitiveVol(TGeoVolume *vol, Int_t px, Int_t py); @@ -85,6 +91,11 @@ public: Int_t nphi=90, Double_t phimin=0., Double_t phimax=360., Double_t rmin=0., Double_t rmax=9999999, Option_t *option=""); virtual void ModifiedPad() const; + virtual void *MakeBox3DBuffer(const TGeoVolume *vol); + virtual void *MakeTube3DBuffer(const TGeoVolume *vol); + virtual void *MakeTubs3DBuffer(const TGeoVolume *vol); + virtual void *MakeSphere3DBuffer(const TGeoVolume *vol); + virtual void *MakePcon3DBuffer(const TGeoVolume *vol); virtual void Paint(Option_t *option=""); void PaintShape(X3DBuffer *buff, Bool_t rangeView, TGeoHMatrix *glmat); void PaintBox(TGeoShape *shape, Option_t *option="", TGeoHMatrix *glmat=0); diff --git a/geompainter/src/TGeoChecker.cxx b/geompainter/src/TGeoChecker.cxx index 9f10a397f9b..9faaf8b83d1 100644 --- a/geompainter/src/TGeoChecker.cxx +++ b/geompainter/src/TGeoChecker.cxx @@ -1,6 +1,6 @@ -// @(#)root/geom:$Name: $:$Id: TGeoChecker.cxx,v 1.21 2003/01/15 18:43:45 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoChecker.cxx,v 1.22 2003/01/31 16:38:23 brun Exp $ // Author: Andrei Gheata 01/11/01 -// TGeoChecker::CheckGeometry() by Mihaela Gheata +// CheckGeometry(), CheckOverlaps() by Mihaela Gheata /************************************************************************* * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. * @@ -35,6 +35,7 @@ #include "TGeoBBox.h" #include "TGeoManager.h" +#include "TGeoPainter.h" #include "TGeoChecker.h" @@ -246,6 +247,134 @@ void TGeoChecker::CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, D delete [] array1; delete [] array2; } + +//----------------------------------------------------------------------------- +void TGeoChecker::CheckOverlaps(const TGeoVolume *vol, Double_t ovlp, Option_t * /*option*/) const +{ +// Check illegal overlaps for volume VOL within a limit OVLP. + if (vol->GetFinder()) return; + Int_t nd = vol->GetNdaughters(); + if (!nd) return; + // first, test if daughters extrude their container + TGeoShape *shapem = vol->GetShape(); + TGeoNode * node; + TGeoMatrix *matrix; + X3DPoints *buff; + Bool_t extrude; + Bool_t ismany; + Double_t *points; + Double_t local[3]; + Double_t point[3]; + Double_t safety = TGeoShape::kBig; + Int_t id, ip; + for (id=0; id<nd; id++) { + node = vol->GetNode(id); + buff = (X3DPoints*)node->GetVolume()->Make3DBuffer(); + if (!buff) { + Error("CheckOverlaps", "could not fill X3D buffer for node %s", node->GetName()); + return; + } + matrix = node->GetMatrix(); + ismany = node->IsOverlapping(); + points = buff->points; + // loop all points + for (ip=0; ip<buff->numPoints; ip++) { + memcpy(local, &points[3*ip], 3*sizeof(Double_t)); + matrix->LocalToMaster(local, point); + extrude = !shapem->Contains(point); + if (extrude) { + safety = shapem->Safety(point, kFALSE); + if (safety<ovlp) extrude=kFALSE; + } + if (extrude) { + if (!ismany) { + printf("=== volume %s : extruding ONLY node %s, extr=%g\n", vol->GetName(), node->GetName(), safety); + } else { + printf("=== volume %s : extruding MANY node %s (G3 intersection), extr=%g\n", vol->GetName(), node->GetName(), safety); + } + break; + } + } + if (points) delete [] points; + delete buff; + } + // now check if the daughters overlap with each other + if (nd<2) return; + TGeoVoxelFinder *vox = vol->GetVoxels(); + TGeoNode *node1; + TGeoMatrix *matrix1; + TGeoShape *shape1; + Bool_t ismany1, overlap; + Int_t novlp; + Int_t *ovlps; + Int_t ko, io; + for (id=0; id<nd; id++) { + node = vol->GetNode(id); + ismany = node->IsOverlapping(); + if (ismany) continue; + shapem = node->GetVolume()->GetShape(); + matrix = node->GetMatrix(); + if (vox) { + vox->FindOverlaps(id); + ovlps = node->GetOverlaps(novlp); + if (!ovlps) continue; + } else continue; + for (ko=0; ko<novlp; ko++) { + io = ovlps[ko]; + if (io<id) continue; + node1 = vol->GetNode(io); + ismany1 = node1->IsOverlapping(); + if (ismany1) continue; + matrix1 = node1->GetMatrix(); + buff = (X3DPoints*)node1->GetVolume()->Make3DBuffer(); + if (!buff) continue; + points = buff->points; + // loop all points + overlap = kFALSE; + for (ip=0; ip<buff->numPoints; ip++) { + memcpy(local, &points[3*ip], 3*sizeof(Double_t)); + matrix1->LocalToMaster(local, point); + matrix->MasterToLocal(point, local); // now point in local reference of node + overlap = shapem->Contains(local); + if (overlap) { + safety = shapem->Safety(local, kTRUE); + if (safety<ovlp) overlap=kFALSE; + } + if (overlap) { + printf("=== volume %s : nodes %s and %s overlapping; overlap=%g\n", + vol->GetName(), node->GetName(), node1->GetName(), safety); + break; + } + } + if (points) delete [] points; + delete buff; + if (overlap) continue; + shape1 = node1->GetVolume()->GetShape(); + buff = (X3DPoints*)node->GetVolume()->Make3DBuffer(); + if (!buff) continue; + points = buff->points; + // loop all points + for (ip=0; ip<buff->numPoints; ip++) { + memcpy(local, &points[3*ip], 3*sizeof(Double_t)); + matrix->LocalToMaster(local, point); + matrix1->MasterToLocal(point, local); // now point in local reference of node + overlap = shape1->Contains(local); + if (overlap) { + safety = shape1->Safety(local, kTRUE); + if (safety<ovlp) overlap=kFALSE; + } + if (overlap) { + printf("=== volume %s : nodes %s and %s overlapping; overlap=%g\n", + vol->GetName(), node->GetName(), node1->GetName(), safety); + break; + } + } + if (points) delete [] points; + delete buff; + } + } +} + //----------------------------------------------------------------------------- void TGeoChecker::CheckPoint(Double_t x, Double_t y, Double_t z, Option_t *) { @@ -269,7 +398,7 @@ void TGeoChecker::CheckPoint(Double_t x, Double_t y, Double_t z, Option_t *) TGeoNode *old = fVsafe->GetNode("SAFETY_1"); if (old) fVsafe->GetNodes()->RemoveAt(vol->GetNdaughters()-1); } - if (vol != fGeom->GetMasterVolume()) fGeom->RestoreMasterVolume(); +// if (vol != fGeom->GetMasterVolume()) fGeom->RestoreMasterVolume(); TGeoNode *node = fGeom->FindNode(point[0], point[1], point[2]); fGeom->MasterToLocal(point, local); Double_t r = TMath::Sqrt(local[0]*local[0]+local[1]*local[1]+local[2]*local[2]); diff --git a/geompainter/src/TGeoPainter.cxx b/geompainter/src/TGeoPainter.cxx index 6ddd5568705..e7fe2780254 100644 --- a/geompainter/src/TGeoPainter.cxx +++ b/geompainter/src/TGeoPainter.cxx @@ -100,6 +100,14 @@ void TGeoPainter::CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, D { fChecker->CheckGeometry(nrays, startx, starty, startz); } + +//______________________________________________________________________________ +void TGeoPainter::CheckOverlaps(const TGeoVolume *vol, Double_t ovlp, Option_t *option) const +{ +// Check overlaps for the top volume of the geometry, within a limit OVLP. + fChecker->CheckOverlaps(vol, ovlp, option); +} + //______________________________________________________________________________ void TGeoPainter::CheckPoint(Double_t x, Double_t y, Double_t z, Option_t *option) { @@ -526,6 +534,25 @@ void TGeoPainter::PaintShape(X3DBuffer *buff, Bool_t rangeView, TGeoHMatrix *glm if (view->GetAutoRange()) view->SetRange(x0,y0,z0,x1,y1,z1,kExpandView); } } + +//______________________________________________________________________________ +void *TGeoPainter::MakeBox3DBuffer(const TGeoVolume *vol) +{ +// Create a box 3D buffer for a given shape. + X3DPoints *buff = new X3DPoints; + const Int_t numpoints = 8; + + buff->numPoints = 8; + + Double_t *points = new Double_t[3*numpoints]; + TGeoShape *shape = vol->GetShape(); + + shape->SetPoints(points); + + buff->points = points; + return buff; +} + //______________________________________________________________________________ void TGeoPainter::PaintBox(TGeoShape *shape, Option_t *option, TGeoHMatrix *glmat) { @@ -621,6 +648,24 @@ void TGeoPainter::PaintCompositeShape(TGeoVolume *vol, Option_t *option) PaintBox(vol->GetShape(), option); } +//______________________________________________________________________________ +void *TGeoPainter::MakeTube3DBuffer(const TGeoVolume *vol) +{ +// Create a box 3D buffer for a given shape. + X3DPoints *buff = new X3DPoints; + Int_t n = fNsegments; + const Int_t numpoints = 4*n; + + buff->numPoints = numpoints; + + Double_t *points = new Double_t[3*numpoints]; + TGeoShape *shape = vol->GetShape(); + + shape->SetPoints(points); + + buff->points = points; + return buff; +} //______________________________________________________________________________ void TGeoPainter::PaintTube(TGeoShape *shape, Option_t *option, TGeoHMatrix *glmat) { @@ -737,6 +782,27 @@ void TGeoPainter::PaintTube(TGeoShape *shape, Option_t *option, TGeoHMatrix *glm if (buff->polys) delete [] buff->polys; if (buff) delete buff; } + +//______________________________________________________________________________ +void *TGeoPainter::MakeTubs3DBuffer(const TGeoVolume *vol) +{ +// Create a box 3D buffer for a given shape. + X3DPoints *buff = new X3DPoints; + + const Int_t n = fNsegments+1; + const Int_t numpoints = 4*n; + + //*-* Allocate memory for points *-* + + Double_t *points = new Double_t[3*numpoints]; + + buff->numPoints = numpoints; + + TGeoShape *shape = vol->GetShape(); + shape->SetPoints(points); + buff->points = points; + return buff; +} //______________________________________________________________________________ void TGeoPainter::PaintTubs(TGeoShape *shape, Option_t *option, TGeoHMatrix *glmat) { @@ -859,11 +925,37 @@ void TGeoPainter::PaintTubs(TGeoShape *shape, Option_t *option, TGeoHMatrix *glm if (buff->polys) delete [] buff->polys; if (buff) delete buff; } +//______________________________________________________________________________ +void *TGeoPainter::MakeSphere3DBuffer(const TGeoVolume *vol) +{ +// Create a box 3D buffer for a given shape. + X3DPoints *buff = new X3DPoints; + + TGeoShape *shape = vol->GetShape(); + ((TGeoSphere*)shape)->SetNumberOfDivisions(fNsegments); + const Int_t n = ((TGeoSphere*)shape)->GetNumberOfDivisions()+1; + Int_t nz = ((TGeoSphere*)shape)->GetNz()+1; + if (nz < 2) return 0; + Int_t numpoints = 2*n*nz; + if (numpoints <= 0) return 0; + + //*-* Allocate memory for points *-* + + Double_t *points = new Double_t[3*numpoints]; + + buff->numPoints = numpoints; + + shape->SetPoints(points); + buff->points = points; + return buff; +} + //______________________________________________________________________________ void TGeoPainter::PaintSphere(TGeoShape *shape, Option_t *option, TGeoHMatrix *glmat) { // paint a sphere Int_t i, j; + ((TGeoSphere*)shape)->SetNumberOfDivisions(fNsegments); const Int_t n = ((TGeoSphere*)shape)->GetNumberOfDivisions()+1; Double_t ph1 = ((TGeoSphere*)shape)->GetPhi1(); Double_t ph2 = ((TGeoSphere*)shape)->GetPhi2(); @@ -1061,6 +1153,26 @@ void TGeoPainter::PaintSphere(TGeoShape *shape, Option_t *option, TGeoHMatrix *g if (buff->polys) delete [] buff->polys; if (buff) delete buff; } + +//______________________________________________________________________________ +void *TGeoPainter::MakePcon3DBuffer(const TGeoVolume *vol) +{ +// Create a box 3D buffer for a given shape. + X3DPoints *buff = new X3DPoints; + + TGeoShape *shape = vol->GetShape(); + const Int_t n = ((TGeoPcon*)shape)->GetNsegments()+1; + Int_t nz = ((TGeoPcon*)shape)->GetNz(); + if (nz < 2) return 0; + Int_t numpoints = nz*2*n; + if (numpoints <= 0) return 0; + Double_t *points = new Double_t[3*numpoints]; + shape->SetPoints(points); + buff->numPoints = numpoints; + buff->points = points; + return buff; +} + //______________________________________________________________________________ void TGeoPainter::PaintPcon(TGeoShape *shape, Option_t *option, TGeoHMatrix *glmat) { -- GitLab