From 4ed4c3bec51a380b511ee239451f3f1a433b7065 Mon Sep 17 00:00:00 2001 From: Rene Brun <Rene.Brun@cern.ch> Date: Wed, 11 Dec 2002 17:10:20 +0000 Subject: [PATCH] From Andrei: attached the last version, with bug fixes in tracking. Now you can do lego plots for everything in alice. It does not do extra steps (except very seldom in ITS in som very thin tubes - 1 extra step to cross). Much more stable. I temporarly disabled TGeoVolume::CheckGeometry from the context menu of volumes since it seems to have a bug that I need to fix. Much better version in computing distances. git-svn-id: http://root.cern.ch/svn/root/trunk@5798 27541ba8-7e3a-0410-8455-c3a389f83636 --- geom/inc/TGeoCache.h | 4 +- geom/inc/TGeoManager.h | 3 +- geom/inc/TGeoVolume.h | 4 +- geom/src/TGeoArb8.cxx | 35 +++++++------- geom/src/TGeoCache.cxx | 17 ++++++- geom/src/TGeoManager.cxx | 7 ++- geom/src/TGeoPgon.cxx | 22 ++++++--- geom/src/TGeoVolume.cxx | 5 +- geompainter/src/TGeoChecker.cxx | 82 +++++++++++---------------------- 9 files changed, 91 insertions(+), 88 deletions(-) diff --git a/geom/inc/TGeoCache.h b/geom/inc/TGeoCache.h index 998f74d0d99..b4bd8c2323a 100644 --- a/geom/inc/TGeoCache.h +++ b/geom/inc/TGeoCache.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoCache.h,v 1.7 2002/12/03 16:01:38 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoCache.h,v 1.8 2002/12/03 17:42:59 rdm Exp $ // Author: Andrei Gheata 18/03/02 /************************************************************************* @@ -155,6 +155,7 @@ public: virtual Int_t GetNfree() const {return (fSize-fNused);} virtual Int_t GetNused() const {return fNused;} virtual const char *GetPath(); + virtual Int_t GetNodeId() const; Int_t GetSize() const {return fSize;} virtual Int_t GetUsageCount() const; virtual void IncreasePool(Int_t size) {fSize+=size;} @@ -219,6 +220,7 @@ public: virtual Int_t GetNfree() const {return fGeoCacheMaxSize;} virtual Int_t GetNused() const {return 0;} virtual const char *GetPath(); + virtual Int_t GetNodeId() const; virtual Int_t GetUsageCount() const {return 0;} virtual void IncreasePool(Int_t /*size*/) {;} virtual void IncrementUsageCount() {;} diff --git a/geom/inc/TGeoManager.h b/geom/inc/TGeoManager.h index 30e6fc8964c..c67dd11386c 100644 --- a/geom/inc/TGeoManager.h +++ b/geom/inc/TGeoManager.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoManager.h,v 1.11 2002/10/21 15:21:13 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoManager.h,v 1.12 2002/10/22 07:43:12 brun Exp $ // Author: Andrei Gheata 25/10/01 /************************************************************************* @@ -270,6 +270,7 @@ public: //--- modeler state getters/setters TGeoNode *GetNode(Int_t level) const {return (TGeoNode*)fNodes->At(level);} + Int_t GetNodeId() const {return fCache->GetNodeId();} TGeoNode *GetNextNode() const {return fNextNode;} TGeoNode *GetMother(Int_t up=1) const {return fCache->GetMother(up);} TGeoHMatrix *GetHMatrix(); diff --git a/geom/inc/TGeoVolume.h b/geom/inc/TGeoVolume.h index 9e03f0696bf..635005aeba1 100644 --- a/geom/inc/TGeoVolume.h +++ b/geom/inc/TGeoVolume.h @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoVolume.h,v 1.12 2002/12/03 16:01:39 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoVolume.h,v 1.13 2002/12/03 17:42:59 rdm Exp $ // Author: Andrei Gheata 30/05/02 /************************************************************************* @@ -94,7 +94,7 @@ public: void ClearNodes() {fNodes = 0;} void ClearShape(); void CleanAll(); - void CheckGeometry(Int_t nrays=1, Double_t startx=0, Double_t starty=0, Double_t startz=0) const; // *MENU* + void CheckGeometry(Int_t nrays=1, Double_t startx=0, Double_t starty=0, Double_t startz=0) const; Int_t CountNodes(Int_t nlevels=1000) const; // *MENU* Bool_t Contains(Double_t *point) const {return fShape->Contains(point);} Bool_t IsFolder() const; diff --git a/geom/src/TGeoArb8.cxx b/geom/src/TGeoArb8.cxx index 2684b89e1b7..091c140d61a 100644 --- a/geom/src/TGeoArb8.cxx +++ b/geom/src/TGeoArb8.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoArb8.cxx,v 1.11 2002/12/03 16:01:39 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoArb8.cxx,v 1.12 2002/12/10 14:34:50 brun Exp $ // Author: Andrei Gheata 31/01/02 /************************************************************************* @@ -300,12 +300,11 @@ Double_t TGeoArb8::DistToPlane(Double_t *point, Double_t *dir, Int_t ipl, Bool_t x1=xs1+tx1*dir[2]*s; x2=xs2+tx2*dir[2]*s; xp=point[0]+s*dir[0]; - if ((x1==x2) || ((xp-x1)*(x2-xp)>=0)) { - y1=ys1+ty1*dir[2]*s; - y2=ys2+ty2*dir[2]*s; - yp=point[1]+s*dir[1]; - if ((y1==y2) || ((yp-y1)*(y2-yp)>=0)) return s; - } + y1=ys1+ty1*dir[2]*s; + y2=ys2+ty2*dir[2]*s; + yp=point[1]+s*dir[1]; + zi = (xp-x1)*(xp-x2)+(yp-y1)*(yp-y2); + if (zi<=0) return s; } } return kBig; @@ -322,12 +321,11 @@ Double_t TGeoArb8::DistToPlane(Double_t *point, Double_t *dir, Int_t ipl, Bool_t x1=xs1+tx1*dir[2]*s; x2=xs2+tx2*dir[2]*s; xp=point[0]+s*dir[0]; - if ((x1==x2) || ((xp-x1)*(x2-xp)>=0)) { - y1=ys1+ty1*dir[2]*s; - y2=ys2+ty2*dir[2]*s; - yp=point[1]+s*dir[1]; - if ((y1==y2) || ((yp-y1)*(y2-yp)>=0)) return s; - } + y1=ys1+ty1*dir[2]*s; + y2=ys2+ty2*dir[2]*s; + yp=point[1]+s*dir[1]; + zi = (xp-x1)*(xp-x2)+(yp-y1)*(yp-y2); + if (zi<=0) return s; } } s=-b+TMath::Sqrt(d); @@ -338,12 +336,11 @@ Double_t TGeoArb8::DistToPlane(Double_t *point, Double_t *dir, Int_t ipl, Bool_t x1=xs1+tx1*dir[2]*s; x2=xs2+tx2*dir[2]*s; xp=point[0]+s*dir[0]; - if ((x1==x2) || ((xp-x1)*(x2-xp)>=0)) { - y1=ys1+ty1*dir[2]*s; - y2=ys2+ty2*dir[2]*s; - yp=point[1]+s*dir[1]; - if ((y1==y2) || ((yp-y1)*(y2-yp)>=0)) return s; - } + y1=ys1+ty1*dir[2]*s; + y2=ys2+ty2*dir[2]*s; + yp=point[1]+s*dir[1]; + zi = (xp-x1)*(xp-x2)+(yp-y1)*(yp-y2); + if (zi<=0) return s; } } } diff --git a/geom/src/TGeoCache.cxx b/geom/src/TGeoCache.cxx index f8aee5b5e95..a2b98cca7d7 100644 --- a/geom/src/TGeoCache.cxx +++ b/geom/src/TGeoCache.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoCache.cxx,v 1.8 2002/10/11 16:41:53 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoCache.cxx,v 1.9 2002/12/03 16:01:39 brun Exp $ // Author: Andrei Gheata 18/03/02 /************************************************************************* @@ -334,6 +334,12 @@ TGeoNode *TGeoNodeCache::GetMother(Int_t up) const return mother; } //----------------------------------------------------------------------------- +Int_t TGeoNodeCache::GetNodeId() const +{ +// Get unique node id. + return fBranch[fLevel]; +} +//----------------------------------------------------------------------------- const char *TGeoNodeCache::GetPath() { // prints the current path @@ -490,6 +496,15 @@ void TGeoCacheDummy::CdUp() fNode = fNodeBranch[fLevel]; fMatrix = fMatrixBranch[fLevel]; } +//----------------------------------------------------------------------------- +Int_t TGeoCacheDummy::GetNodeId() const +{ +// Get unique node id. + Int_t id=0; + for (Int_t level=0;level<fLevel+1; level++) id += (Int_t)fNodeBranch[level]; + return id; +} + //----------------------------------------------------------------------------- const char *TGeoCacheDummy::GetPath() { diff --git a/geom/src/TGeoManager.cxx b/geom/src/TGeoManager.cxx index 49224ad3370..c8c570fcc4f 100644 --- a/geom/src/TGeoManager.cxx +++ b/geom/src/TGeoManager.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoManager.cxx,v 1.23 2002/11/15 14:41:19 rdm Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoManager.cxx,v 1.24 2002/12/03 16:01:39 brun Exp $ // Author: Andrei Gheata 25/10/01 /************************************************************************* @@ -2147,6 +2147,7 @@ TGeoNode *TGeoManager::Step(Bool_t is_geom, Bool_t cross) } if (is_geom) epsil=(cross)?1E-6:-1E-6; TGeoNode *old = fCurrentNode; + Int_t idold = GetNodeId(); if (fIsOutside) old = 0; fStep += epsil; for (Int_t i=0; i<3; i++) fPoint[i]+=fStep*fDirection[i]; @@ -2154,6 +2155,10 @@ TGeoNode *TGeoManager::Step(Bool_t is_geom, Bool_t cross) if (fIsOutside) current=0; if (is_geom) { fIsEntering = (current==old)?kFALSE:kTRUE; + if (!fIsEntering) { + Int_t id = GetNodeId(); + fIsEntering = (id==idold)?kFALSE:kTRUE; + } fIsExiting = !fIsEntering; if (fIsEntering && fIsNullStep) fIsNullStep = kFALSE; fIsOnBoundary = kTRUE; diff --git a/geom/src/TGeoPgon.cxx b/geom/src/TGeoPgon.cxx index f72fc86aa33..ff549f52f89 100644 --- a/geom/src/TGeoPgon.cxx +++ b/geom/src/TGeoPgon.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoPgon.cxx,v 1.9 2002/12/06 16:45:03 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoPgon.cxx,v 1.10 2002/12/10 14:34:50 brun Exp $ // Author: Andrei Gheata 31/01/02 // TGeoPgon::Contains() implemented by Mihaela Gheata @@ -453,12 +453,20 @@ Double_t TGeoPgon::DistToOutSect(Double_t *point, Double_t *dir, Int_t &iz, Int_ { // compute distance to outside from a pgon phi trapezoid Double_t saf; + Double_t dmin = kBig; Double_t snext[6]; Int_t i; for (i=0; i<6; i++) snext[i]=kBig; Double_t zmin = fZ[iz]; Double_t zmax = fZ[iz+1]; - if (zmax==zmin) return 1E-12; + if (zmax==zmin) { + iz += (dir[2]>0)?1:-1; + Double_t pt[3]; + for (i=0; i<3; i++) pt[i] = point[i]+1E-8*dir[i]; + dmin = 0.; + if (Contains(&pt[0])) dmin = DistToOutSect(&pt[0], dir, iz, isect); + return (dmin+1E-8); + } Double_t divphi = fDphi/fNedges; Double_t phi1 = (fPhi1 + divphi*(isect-1))*kDegRad; Double_t phi2 = phi1 + divphi*kDegRad; @@ -466,7 +474,6 @@ Double_t TGeoPgon::DistToOutSect(Double_t *point, Double_t *dir, Int_t &iz, Int_ Double_t cphim = TMath::Cos(phim); Double_t sphim = TMath::Sin(phim); Double_t minsafe = 0; - Double_t dmin = kBig; Double_t no[3]; // check outer slanted face Double_t ct, st; @@ -534,8 +541,9 @@ Double_t TGeoPgon::DistToOutSect(Double_t *point, Double_t *dir, Int_t &iz, Int_ // z plane crossed iz += 2*icheck-5; if ((iz<0) || (iz>(fNz-2))) return dmin; - for (i=0; i<3; i++) pt[i]=point[i]+(dmin+1E-10)*dir[i]; - dmin += DistToOutSect(&pt[0], dir, iz, isect)+1E-10; + for (i=0; i<3; i++) pt[i]=point[i]+(dmin+1E-8)*dir[i]; + if (Contains(&pt[0])) + dmin += DistToOutSect(&pt[0], dir, iz, isect)+1E-8; return dmin; } isect += 2*icheck-9; @@ -545,8 +553,8 @@ Double_t TGeoPgon::DistToOutSect(Double_t *point, Double_t *dir, Int_t &iz, Int_ } else { if ((isect<1) || (isect>fNedges)) return dmin; } - for (i=0; i<3; i++) pt[i]=point[i]+(dmin+1E-10)*dir[i]; - dmin += DistToOutSect(&pt[0], dir, iz, isect)+1E-10; + for (i=0; i<3; i++) pt[i]=point[i]+(dmin+1E-8)*dir[i]; + dmin += DistToOutSect(&pt[0], dir, iz, isect)+1E-8; return dmin; } //----------------------------------------------------------------------------- diff --git a/geom/src/TGeoVolume.cxx b/geom/src/TGeoVolume.cxx index 010785d07f1..b6ffb18fa8c 100644 --- a/geom/src/TGeoVolume.cxx +++ b/geom/src/TGeoVolume.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoVolume.cxx,v 1.13 2002/12/03 16:01:40 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoVolume.cxx,v 1.14 2002/12/10 07:52:33 brun Exp $ // Author: Andrei Gheata 30/05/02 // Divide() implemented by Mihaela Gheata @@ -410,6 +410,9 @@ void TGeoVolume::Draw(Option_t *option) void TGeoVolume::DrawOnly(Option_t *option) { // draw only this volume + TGeoVolume *old_vol = gGeoManager->GetTopVolume(); + if (old_vol!=this) gGeoManager->SetTopVolume(this); + else old_vol=0; TVirtualGeoPainter *painter = gGeoManager->GetGeomPainter(); if (!painter) return; painter->DrawOnly(option); diff --git a/geompainter/src/TGeoChecker.cxx b/geompainter/src/TGeoChecker.cxx index 7af6deeba1a..d408b2b5e9f 100644 --- a/geompainter/src/TGeoChecker.cxx +++ b/geompainter/src/TGeoChecker.cxx @@ -1,4 +1,4 @@ -// @(#)root/geom:$Name: $:$Id: TGeoChecker.cxx,v 1.17 2002/12/10 08:40:15 brun Exp $ +// @(#)root/geom:$Name: $:$Id: TGeoChecker.cxx,v 1.18 2002/12/10 14:34:50 brun Exp $ // Author: Andrei Gheata 01/11/01 // TGeoChecker::CheckGeometry() by Mihaela Gheata @@ -79,8 +79,8 @@ void TGeoChecker::CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, D Double_t dir[3]; Double_t dummy[3]; Double_t eps = 0.; - Double_t *array1 = new Double_t[3*100]; - Double_t *array2 = new Double_t[3*100]; + Double_t *array1 = new Double_t[3*1000]; + Double_t *array2 = new Double_t[3*1000]; TObjArray *pma = new TObjArray(); TPolyMarker3D *pm; pm = new TPolyMarker3D(); @@ -99,7 +99,7 @@ void TGeoChecker::CheckGeometry(Int_t nrays, Double_t startx, Double_t starty, D pm->SetMarkerSize(0.4); pma->AddAt(pm, 2); Int_t nelem1, nelem2; - Int_t dim1=100, dim2=100; + Int_t dim1=1000, dim2=1000; if ((startx==0) && (starty==0) && (startz==0)) eps=1E-3; start[0] = startx+eps; start[1] = starty+eps; @@ -284,8 +284,10 @@ void TGeoChecker::CheckPoint(Double_t x, Double_t y, Double_t z, Option_t *) pm->SetNextPoint(local[0], local[1], local[2]); if (vol->GetNdaughters()) { fGeom->SetVisLevel(1); + if (!vol->IsVisible()) vol->SetVisibility(kTRUE); vol->Draw(); } else { + if (!vol->IsVisible()) vol->SetVisibility(kTRUE); vol->DrawOnly(); } pm->Draw("SAME"); @@ -459,14 +461,13 @@ void TGeoChecker::RandomRays(Int_t nrays, Double_t startx, Double_t starty, Doub Double_t start[3]; Double_t dir[3]; Double_t eps = 0.; + Int_t istep= 0; if ((startx==0) && (starty==0) && (startz==0)) eps=1E-3; Double_t *point = fGeom->GetCurrentPoint(); vol->Draw(); printf("Start... %i rays\n", nrays); - //TGeoNode *node; - //Bool_t is_sentering; TGeoNode *startnode, *endnode; - Bool_t vis1,vis2, is_entering; + Bool_t vis1,vis2; Int_t i=0; Int_t ipoint; Int_t itot=0; @@ -481,7 +482,6 @@ void TGeoChecker::RandomRays(Int_t nrays, Double_t startx, Double_t starty, Doub start[0] = startx+eps; start[1] = starty+eps; start[2] = startz+eps; -// printf("startpoint : %g, %g, %g\n", startx, starty, startz); phi = 2*TMath::Pi()*gRandom->Rndm(); theta= TMath::ACos(1.-2.*gRandom->Rndm()); dir[0]=TMath::Sin(theta)*TMath::Cos(phi); @@ -490,82 +490,54 @@ void TGeoChecker::RandomRays(Int_t nrays, Double_t startx, Double_t starty, Doub fGeom->InitTrack(&start[0], &dir[0]); line = 0; startnode = fGeom->GetCurrentNode(); -// if (startnode) printf("---first : %s\n", startnode->GetName()); -// else printf("---first : NULL\n"); if (fGeom->IsOutside()) startnode=0; vis1 = (startnode)?(startnode->IsOnScreen()):kFALSE; if (vis1) { -// printf(" >>>starting first segment\n"); -// if (startnode) printf( ">>>start = %s\n", startnode->GetName()); line = new TPolyLine3D(2); line->SetLineColor(startnode->GetVolume()->GetLineColor()); line->SetPoint(ipoint++, startx, starty, startz); i++; pm->Add(line); - } //else printf(" invisible, no segment\n"); - // find the node that will be crossed first - fGeom->FindNextBoundary(); - fGeom->IsStepEntering(); + } // find where we end-up - endnode = fGeom->Step(); - if (fGeom->IsOutside()) endnode=0; -// if (endnode) printf("---end : %s\n", endnode->GetName()); -// else printf("---end : NULL\n"); + fGeom->FindNextBoundary(); step = fGeom->GetStep(); + endnode = fGeom->Step(); vis2 = (endnode)?(endnode->IsOnScreen()):kFALSE; -// if (endnode) printf(">>>endnode : %s\n", endnode->GetName()); -// if (vis2) printf(" end visible\n"); -// else printf(" end invisible\n"); - is_entering = fGeom->IsEntering(); -// is_null = fGeom->IsNullStep(); -// printf("endpoint : %g, %g, %g step: %g\n", fGeom->GetCurrentPoint()[0], -// fGeom->GetCurrentPoint()[1], fGeom->GetCurrentPoint()[2], step); -// printf("propagating...\n"); while (step<1E10) { - if (is_entering && ipoint>0) { -// printf(" ending segment\n"); + while (!fGeom->IsEntering()) { + istep++; + if (istep>1E4) break; + fGeom->SetStep(1E-3); + endnode = fGeom->Step(); + step += 1E-3; + } + if (istep>1E4) break; +// if (istep) printf("ADDED : %f (%i steps)\n", istep*1E-3, istep); + vis2 = (endnode)?(endnode->IsOnScreen()):kFALSE; + if (ipoint>0) { // old visible node had an entry point -> finish segment line->SetPoint(ipoint, point[0], point[1], point[2]); ipoint = 0; line = 0; } - if (is_entering && vis2) { + if (vis2) { // create new segment -// printf(" >>>new segment for : %s\n", endnode->GetName()); line = new TPolyLine3D(2); line->SetLineColor(endnode->GetVolume()->GetLineColor()); line->SetPoint(ipoint++, point[0], point[1], point[2]); i++; pm->Add(line); - } // else printf(" entering=%i vis2=%i - no segment\n", (Int_t)is_entering, (Int_t)vis2); + } // now see if we can make an other step - if (endnode==0 && step>1E10) { -// printf("NULL. End track.\n"); - break; - } - while (!fGeom->IsEntering()) { - fGeom->SetStep(1E-3); - fGeom->Step(); - } -// if (is_null) printf("null step, start:%f %f %f\n", fGeom->GetCurrentPoint()[0], -// fGeom->GetCurrentPoint()[1], fGeom->GetCurrentPoint()[2]); + if (endnode==0 && step>1E10) break; + istep = 0; // generate an extra step to cross boundary startnode = endnode; fGeom->FindNextBoundary(); -// fGeom->IsStepEntering(); - endnode = fGeom->Step(); - if (fGeom->IsOutside()) endnode=0; -// if (endnode) printf(">>>new node: %s\n", endnode->GetName()); step = fGeom->GetStep(); + endnode = fGeom->Step(); vis2 = (endnode)?(endnode->IsOnScreen()):kFALSE; - is_entering = fGeom->IsEntering(); -// is_null = fGeom->IsNullStep(); -// printf(" new step : %g point: %g, %g, %g\n", step, fGeom->GetCurrentPoint()[0], -// fGeom->GetCurrentPoint()[1], fGeom->GetCurrentPoint()[2]); -// if (is_entering) printf(" entering\n"); -// else printf(" not entering, same node\n"); -// if (endnode) printf(" node after step : %s\n", endnode->GetName()); -// else printf(" NULL node after step\n"); } } // draw all segments -- GitLab