ROOT logo
// @(#)root/geom:$Id: TGeoEltu.cxx 27922 2009-03-24 18:20:59Z brun $
// Author: Mihaela Gheata   05/06/02

/*************************************************************************
 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

//_____________________________________________________________________________
// TGeoEltu - elliptical tube class. It takes 3 parameters : 
// semi-axis of the ellipse along x, semi-asix of the ellipse along y
// and half-length dz. 
//
//_____________________________________________________________________________
//Begin_Html
/*
<img src="gif/t_eltu.gif">
*/
//End_Html

#include "Riostream.h"

#include "TGeoManager.h"
#include "TGeoVolume.h"
#include "TGeoEltu.h"
#include "TBuffer3D.h"
#include "TBuffer3DTypes.h"
#include "TMath.h"

ClassImp(TGeoEltu)
   
//_____________________________________________________________________________
TGeoEltu::TGeoEltu()
{
// Dummy constructor
   SetShapeBit(TGeoShape::kGeoEltu);
}   

//_____________________________________________________________________________
TGeoEltu::TGeoEltu(Double_t a, Double_t b, Double_t dz)
           :TGeoTube()
{
// Default constructor specifying X and Y semiaxis length
   SetShapeBit(TGeoShape::kGeoEltu);
   SetEltuDimensions(a, b, dz);
   ComputeBBox();
}

//_____________________________________________________________________________
TGeoEltu::TGeoEltu(const char *name, Double_t a, Double_t b, Double_t dz)
           :TGeoTube(name,0.,b,dz)
{
// Default constructor specifying X and Y semiaxis length
   SetName(name);
   SetShapeBit(TGeoShape::kGeoEltu);
   SetEltuDimensions(a, b, dz);
   ComputeBBox();
}

//_____________________________________________________________________________
TGeoEltu::TGeoEltu(Double_t *param)
{
// Default constructor specifying minimum and maximum radius
// param[0] =  A
// param[1] =  B
// param[2] = dz
   SetShapeBit(TGeoShape::kGeoEltu);
   SetDimensions(param);
   ComputeBBox();
}

//_____________________________________________________________________________
TGeoEltu::~TGeoEltu()
{
// destructor
}

//_____________________________________________________________________________
Double_t TGeoEltu::Capacity() const
{
// Computes capacity of the shape in [length^3]
   Double_t capacity = 2.*TMath::Pi()*fDz*fRmin*fRmax;
   return capacity;
}   

//_____________________________________________________________________________   
void TGeoEltu::ComputeBBox()
{
// compute bounding box of the tube
   fDX = fRmin;
   fDY = fRmax;
   fDZ = fDz;
}   

//_____________________________________________________________________________   
void TGeoEltu::ComputeNormal(Double_t *point, Double_t *dir, Double_t *norm)
{
// Compute normal to closest surface from POINT.
   Double_t a = fRmin;
   Double_t b = fRmax;
   Double_t eps = TMath::Sqrt(point[0]*point[0]/(a*a)+point[1]*point[1]/(b*b))-1.;
   if (eps<1E-4 && TMath::Abs(fDz-TMath::Abs(point[2]))<1E-5) {
      norm[0] = norm[1] = 0;
      norm[2] = TMath::Sign(1.,dir[2]);
      return;
   }   
   norm[2] = 0.;
   Double_t r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
   Double_t st = point[1]/r;
   Double_t ct = point[0]/r;
   Double_t rr = TMath::Sqrt(b*b*ct*ct+a*a*st*st);
   norm[0] = b*ct/rr;
   norm[1] = a*st/rr;
   if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
      norm[0] = -norm[0];
      norm[1] = -norm[1];
   }   
}

//_____________________________________________________________________________
Bool_t TGeoEltu::Contains(Double_t *point) const
{
// test if point is inside the elliptical tube
   if (TMath::Abs(point[2]) > fDz) return kFALSE;
   Double_t r2 = (point[0]*point[0])/(fRmin*fRmin)+(point[1]*point[1])/(fRmax*fRmax);
   if (r2>1.)  return kFALSE;
   return kTRUE;
}

//_____________________________________________________________________________
Int_t TGeoEltu::DistancetoPrimitive(Int_t px, Int_t py)
{
// compute closest distance from point px,py to each vertex
   Int_t n = gGeoManager->GetNsegments();
   const Int_t numPoints=4*n;
   return ShapeDistancetoPrimitive(numPoints, px, py);
}   

//_____________________________________________________________________________
Double_t TGeoEltu::DistFromInside(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 tube
   Double_t a2=fRmin*fRmin;
   Double_t b2=fRmax*fRmax;
   Double_t safz1=fDz-point[2];
   Double_t safz2=fDz+point[2];
   
   if (iact<3 && safe) {
      Double_t x0=TMath::Abs(point[0]);
      Double_t y0=TMath::Abs(point[1]);
      Double_t x1=x0;
      Double_t y1=TMath::Sqrt((fRmin-x0)*(fRmin+x0))*fRmax/fRmin;
      Double_t y2=y0;
      Double_t x2=TMath::Sqrt((fRmax-y0)*(fRmax+y0))*fRmin/fRmax;
      Double_t d1=(x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
      Double_t d2=(x2-x0)*(x2-x0)+(y2-y0)*(y2-y0);
      Double_t x3,y3;
   
      Double_t safr=TGeoShape::Big();
      Double_t safz = TMath::Min(safz1,safz2);
      for (Int_t i=0; i<8; i++) {
         if (fRmax<fRmin) {
            x3=0.5*(x1+x2);
            y3=TMath::Sqrt((fRmin-x3)*(fRmin+x3))*fRmax/fRmin;;
         } else {
            y3=0.5*(y1+y2);   
            x3=TMath::Sqrt((fRmax-y3)*(fRmax+y3))*fRmin/fRmax;
         }
         if (d1<d2) {
            x2=x3;
            y2=y3;
            d2=(x2-x0)*(x2-x0)+(y2-y0)*(y2-y0);
         } else {
            x1=x3;
            y1=y3;
            d1=(x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
         }
      }
      safr=TMath::Sqrt(d1)-1.0E-3;   
      *safe = TMath::Min(safz, safr);
      if (iact==0) return TGeoShape::Big();
      if ((iact==1) && (*safe>step)) return TGeoShape::Big();
   }
   // compute distance to surface 
   // Do Z
   Double_t snxt = TGeoShape::Big();
   if (dir[2]>0) {
      snxt=safz1/dir[2];
   } else {
      if (dir[2]<0) snxt=-safz2/dir[2];
   } 
   // do eliptical surface
   Double_t sz = snxt;
   Double_t xz=point[0]+dir[0]*sz;
   Double_t yz=point[1]+dir[1]*sz;
   if ((xz*xz/a2+yz*yz/b2)<=1) return snxt;
   Double_t u=dir[0]*dir[0]*b2+dir[1]*dir[1]*a2;
   Double_t v=point[0]*dir[0]*b2+point[1]*dir[1]*a2;
   Double_t w=point[0]*point[0]*b2+point[1]*point[1]*a2-a2*b2;
   Double_t d=v*v-u*w;
   if (d<0) return snxt;
   if (TGeoShape::IsSameWithinTolerance(u,0)) return snxt;
   Double_t sd=TMath::Sqrt(d);
   Double_t tau1=(-v+sd)/u;
   Double_t tau2=(-v+sd)/u;
   
   if (tau1<0) {
      if (tau2<0) return snxt;
   } else {
      snxt=tau1;
      if ((tau2>0) && (tau2<tau1)) return tau2;
   }
   return snxt;      
}

//_____________________________________________________________________________
Double_t TGeoEltu::DistFromOutside(Double_t *point, Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
{
// compute distance from outside point to surface of the tube and safe distance
   Double_t snxt=TGeoShape::Big();
   Double_t safz=TMath::Abs(point[2])-fDz;
   Double_t a2=fRmin*fRmin;
   Double_t b2=fRmax*fRmax;
   if (iact<3 && safe) {
      Double_t x0=TMath::Abs(point[0]);
      Double_t y0=TMath::Abs(point[1]);
      *safe=0.;
      if ((x0*x0/a2+y0*y0/b2)>=1) {
         Double_t phi1=0;
         Double_t phi2=0.5*TMath::Pi();
         Double_t phi3;
         Double_t x3,y3,d;
         for (Int_t i=0; i<10; i++) {
            phi3=(phi1+phi2)*0.5;
            x3=fRmin*TMath::Cos(phi3);
            y3=fRmax*TMath::Sin(phi3);
            d=y3*a2*(x0-x3)-x3*b2*(y0-y3);
            if (d<0) phi1=phi3;
            else phi2=phi3;
         }
         *safe=TMath::Sqrt((x0-x3)*(x0-x3)+(y0-y3)*(y0-y3));
      }
      if (safz>0) {
         *safe=TMath::Sqrt((*safe)*(*safe)+safz*safz);
      } 
      if (iact==0) return TGeoShape::Big();
      if ((iact==1) && (step<*safe)) return TGeoShape::Big();
   }
   // compute vector distance
   if ((safz>0) && (point[2]*dir[2]>=0)) return TGeoShape::Big();
// Check if the bounding box is crossed within the requested distance
   Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
   if (sdist>=step) return TGeoShape::Big();
   Double_t zi;
   if (!TGeoShape::IsSameWithinTolerance(dir[2],0)) {
      Double_t u=dir[0]*dir[0]*b2+dir[1]*dir[1]*a2;
      if (TGeoShape::IsSameWithinTolerance(u,0)) return TGeoShape::Big();
      Double_t v=point[0]*dir[0]*b2+point[1]*dir[1]*a2;
      Double_t w=point[0]*point[0]*b2+point[1]*point[1]*a2-a2*b2;
      Double_t d=v*v-u*w;
      if (d<0) return TGeoShape::Big();
      Double_t dsq=TMath::Sqrt(d);
      Double_t tau[2];
      tau[0]=(-v+dsq)/u;
      tau[1]=(-v-dsq)/u;
      for (Int_t j=0; j<2; j++) {
         if (tau[j]>=0) {
            zi=point[2]+tau[j]*dir[2];
            if ((TMath::Abs(zi)-fDz)<0)
               snxt=TMath::Min(snxt,tau[j]);
         } 
      }
   }   
   // do z
   zi=TGeoShape::Big();
   if (safz>0) {
      if (TGeoShape::IsSameWithinTolerance(dir[2],0)) return TGeoShape::Big();
      if (point[2]>0) zi=fDz;
      if (point[2]<0) zi=-fDz;     
      Double_t tauz=(zi-point[2])/dir[2];
      Double_t xz=point[0]+dir[0]*tauz;
      Double_t yz=point[1]+dir[1]*tauz;
      if ((xz*xz/a2+yz*yz/b2)<=1) snxt=tauz;
   }
   return snxt;   
}

//_____________________________________________________________________________
TGeoVolume *TGeoEltu::Divide(TGeoVolume * /*voldiv*/, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/, 
                             Double_t /*start*/, Double_t /*step*/) 
{
// Divide the shape along one axis.
   Error("Divide", "Elliptical tubes divisions not implemenetd");
   return 0;
}   

//_____________________________________________________________________________
void TGeoEltu::GetBoundingCylinder(Double_t *param) const
{
//--- Fill vector param[4] with the bounding cylinder parameters. The order
// is the following : Rmin, Rmax, Phi1, Phi2
   param[0] = 0.;                  // Rmin
   param[1] = TMath::Max(fRmin, fRmax); // Rmax
   param[1] *= param[1];
   param[2] = 0.;                  // Phi1
   param[3] = 360.;                // Phi2 
}   

//_____________________________________________________________________________
TGeoShape *TGeoEltu::GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix * /*mat*/) const
{
// in case shape has some negative parameters, these has to be computed
// in order to fit the mother
   if (!TestShapeBit(kGeoRunTimeShape)) return 0;
   if (!mother->TestShapeBit(kGeoEltu)) {
      Error("GetMakeRuntimeShape", "invalid mother");
      return 0;
   }
   Double_t a, b, dz;
   a = fRmin;
   b = fRmax;
   dz = fDz;
   if (fDz<0) dz=((TGeoEltu*)mother)->GetDz();
   if (fRmin<0)
      a = ((TGeoEltu*)mother)->GetA();
   if (fRmax<0) 
      a = ((TGeoEltu*)mother)->GetB();

   return (new TGeoEltu(a, b, dz));
}

//_____________________________________________________________________________
void TGeoEltu::InspectShape() const
{
// print shape parameters
   printf("*** Shape %s: TGeoEltu ***\n", GetName());
   printf("    A    = %11.5f\n", fRmin);
   printf("    B    = %11.5f\n", fRmax);
   printf("    dz   = %11.5f\n", fDz);
   printf(" Bounding box:\n");
   TGeoBBox::InspectShape();
}

//_____________________________________________________________________________
Double_t TGeoEltu::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 x0 = TMath::Abs(point[0]);
   Double_t y0 = TMath::Abs(point[1]);
   Double_t x1, y1, dx, dy;
   Double_t safr, safz;
   safr = safz = TGeoShape::Big();
   if (in) {
      x1 = fRmin*TMath::Sqrt(1.-(y0*y0)/(fRmax*fRmax));
      y1 = fRmax*TMath::Sqrt(1.-(x0*x0)/(fRmin*fRmin));
      dx = x1-x0;
      dy = y1-y0;
      if (TMath::Abs(dx)<TGeoShape::Tolerance()) return 0;
      safr = dx*dy/TMath::Sqrt(dx*dx+dy*dy);
      safz = fDz - TMath::Abs(point[2]);
      return TMath::Min(safr,safz);
   }   

   if (TMath::Abs(x0)<TGeoShape::Tolerance()) {
      safr = y0 - fRmax;
   } else {
      if (TMath::Abs(y0)<TGeoShape::Tolerance()) {
         safr = x0 - fRmin;
      } else {
         Double_t f = fRmin*fRmax/TMath::Sqrt(x0*x0*fRmax*fRmax+y0*y0*fRmin*fRmin);
         x1 = f*x0;
         y1 = f*y0;
         dx = x0-x1;
         dy = y0-y1;
         Double_t ast = fRmin*y1/fRmax;
         Double_t bct = fRmax*x1/fRmin;
         Double_t d = TMath::Sqrt(bct*bct+ast*ast);
         safr = (dx*bct+dy*ast)/d;
      }
   }
   safz = TMath::Abs(point[2])-fDz;            
   return TMath::Max(safr, safz);
}

//_____________________________________________________________________________
void TGeoEltu::SavePrimitive(ostream &out, Option_t * /*option*/ /*= ""*/)
{
// Save a primitive as a C++ statement(s) on output stream "out".
   if (TObject::TestBit(kGeoSavePrimitive)) return;
   out << "   // Shape: " << GetName() << " type: " << ClassName() << endl;
   out << "   a  = " << fRmin << ";" << endl;
   out << "   b  = " << fRmax << ";" << endl;
   out << "   dz = " << fDz << ";" << endl;
   out << "   TGeoShape *" << GetPointerName() << " = new TGeoEltu(\"" << GetName() << "\",a,b,dz);" << endl;
   TObject::SetBit(TGeoShape::kGeoSavePrimitive);
}

//_____________________________________________________________________________
void TGeoEltu::SetEltuDimensions(Double_t a, Double_t b, Double_t dz)
{
// Set dimensions of the eliptical tube.
   if ((a<=0) || (b<0) || (dz<0)) {
      SetShapeBit(kGeoRunTimeShape);
   }
   fRmin=a;
   fRmax=b;
   fDz=dz;
}   

//_____________________________________________________________________________
void TGeoEltu::SetDimensions(Double_t *param)
{
// Set shape dimensions starting from an array.
   Double_t a    = param[0];
   Double_t b    = param[1];
   Double_t dz   = param[2];
   SetEltuDimensions(a, b, dz);
}   

//_____________________________________________________________________________
void TGeoEltu::SetPoints(Double_t *points) const
{
// Create eliptical tube mesh points
   Double_t dz;
   Int_t j, n;

   n = gGeoManager->GetNsegments();
   Double_t dphi = 360./n;
   Double_t phi = 0;
   Double_t cph,sph;
   dz = fDz;

   Int_t indx = 0;
   Double_t r2,r;
   Double_t a2=fRmin*fRmin;
   Double_t b2=fRmax*fRmax;

   if (points) {
      for (j = 0; j < n; j++) {
         points[indx+6*n] = points[indx] = 0;
         indx++;
         points[indx+6*n] = points[indx] = 0;
         indx++;
         points[indx+6*n] = dz;
         points[indx]     =-dz;
         indx++;
      }
      for (j = 0; j < n; j++) {
         phi = j*dphi*TMath::DegToRad();
         sph=TMath::Sin(phi);
         cph=TMath::Cos(phi);
         r2=(a2*b2)/(b2+(a2-b2)*sph*sph);
         r=TMath::Sqrt(r2);
         points[indx+6*n] = points[indx] = r*cph;
         indx++;
         points[indx+6*n] = points[indx] = r*sph;
         indx++;
         points[indx+6*n]= dz;
         points[indx]    =-dz;
         indx++;
      }
   }
}

//_____________________________________________________________________________
void TGeoEltu::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
{
// Returns numbers of vertices, segments and polygons composing the shape mesh.
   TGeoTube::GetMeshNumbers(nvert,nsegs,npols);
}

//_____________________________________________________________________________
Int_t TGeoEltu::GetNmeshVertices() const
{
// Returns the number of vertices on the mesh.
   return TGeoTube::GetNmeshVertices();
}   
   
//_____________________________________________________________________________
void TGeoEltu::SetPoints(Float_t *points) const
{
// Create eliptical tube mesh points
   Double_t dz;
   Int_t j, n;

   n = gGeoManager->GetNsegments();
   Double_t dphi = 360./n;
   Double_t phi = 0;
   Double_t cph,sph;
   dz = fDz;

   Int_t indx = 0;
   Double_t r2,r;
   Double_t a2=fRmin*fRmin;
   Double_t b2=fRmax*fRmax;

   if (points) {
      for (j = 0; j < n; j++) {
         points[indx+6*n] = points[indx] = 0;
         indx++;
         points[indx+6*n] = points[indx] = 0;
         indx++;
         points[indx+6*n] = dz;
         points[indx]     =-dz;
         indx++;
      }
      for (j = 0; j < n; j++) {
         phi = j*dphi*TMath::DegToRad();
         sph=TMath::Sin(phi);
         cph=TMath::Cos(phi);
         r2=(a2*b2)/(b2+(a2-b2)*sph*sph);
         r=TMath::Sqrt(r2);
         points[indx+6*n] = points[indx] = r*cph;
         indx++;
         points[indx+6*n] = points[indx] = r*sph;
         indx++;
         points[indx+6*n]= dz;
         points[indx]    =-dz;
         indx++;
      }
   }
}

//_____________________________________________________________________________
const TBuffer3D & TGeoEltu::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
{
// Fills a static 3D buffer and returns a reference.
   static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
   TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);

   if (reqSections & TBuffer3D::kRawSizes) {
      Int_t n = gGeoManager->GetNsegments();
      Int_t nbPnts = 4*n;
      Int_t nbSegs = 8*n;
      Int_t nbPols = 4*n;      
      if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
         buffer.SetSectionsValid(TBuffer3D::kRawSizes);
      }
   }
   if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
      SetPoints(buffer.fPnts);
      if (!buffer.fLocalFrame) {
         TransformPoints(buffer.fPnts, buffer.NbPnts());
      }
      SetSegsAndPols(buffer);  
      buffer.SetSectionsValid(TBuffer3D::kRaw);
   }
      
   return buffer;
}
 TGeoEltu.cxx:1
 TGeoEltu.cxx:2
 TGeoEltu.cxx:3
 TGeoEltu.cxx:4
 TGeoEltu.cxx:5
 TGeoEltu.cxx:6
 TGeoEltu.cxx:7
 TGeoEltu.cxx:8
 TGeoEltu.cxx:9
 TGeoEltu.cxx:10
 TGeoEltu.cxx:11
 TGeoEltu.cxx:12
 TGeoEltu.cxx:13
 TGeoEltu.cxx:14
 TGeoEltu.cxx:15
 TGeoEltu.cxx:16
 TGeoEltu.cxx:17
 TGeoEltu.cxx:18
 TGeoEltu.cxx:19
 TGeoEltu.cxx:20
 TGeoEltu.cxx:21
 TGeoEltu.cxx:22
 TGeoEltu.cxx:23
 TGeoEltu.cxx:24
 TGeoEltu.cxx:25
 TGeoEltu.cxx:26
 TGeoEltu.cxx:27
 TGeoEltu.cxx:28
 TGeoEltu.cxx:29
 TGeoEltu.cxx:30
 TGeoEltu.cxx:31
 TGeoEltu.cxx:32
 TGeoEltu.cxx:33
 TGeoEltu.cxx:34
 TGeoEltu.cxx:35
 TGeoEltu.cxx:36
 TGeoEltu.cxx:37
 TGeoEltu.cxx:38
 TGeoEltu.cxx:39
 TGeoEltu.cxx:40
 TGeoEltu.cxx:41
 TGeoEltu.cxx:42
 TGeoEltu.cxx:43
 TGeoEltu.cxx:44
 TGeoEltu.cxx:45
 TGeoEltu.cxx:46
 TGeoEltu.cxx:47
 TGeoEltu.cxx:48
 TGeoEltu.cxx:49
 TGeoEltu.cxx:50
 TGeoEltu.cxx:51
 TGeoEltu.cxx:52
 TGeoEltu.cxx:53
 TGeoEltu.cxx:54
 TGeoEltu.cxx:55
 TGeoEltu.cxx:56
 TGeoEltu.cxx:57
 TGeoEltu.cxx:58
 TGeoEltu.cxx:59
 TGeoEltu.cxx:60
 TGeoEltu.cxx:61
 TGeoEltu.cxx:62
 TGeoEltu.cxx:63
 TGeoEltu.cxx:64
 TGeoEltu.cxx:65
 TGeoEltu.cxx:66
 TGeoEltu.cxx:67
 TGeoEltu.cxx:68
 TGeoEltu.cxx:69
 TGeoEltu.cxx:70
 TGeoEltu.cxx:71
 TGeoEltu.cxx:72
 TGeoEltu.cxx:73
 TGeoEltu.cxx:74
 TGeoEltu.cxx:75
 TGeoEltu.cxx:76
 TGeoEltu.cxx:77
 TGeoEltu.cxx:78
 TGeoEltu.cxx:79
 TGeoEltu.cxx:80
 TGeoEltu.cxx:81
 TGeoEltu.cxx:82
 TGeoEltu.cxx:83
 TGeoEltu.cxx:84
 TGeoEltu.cxx:85
 TGeoEltu.cxx:86
 TGeoEltu.cxx:87
 TGeoEltu.cxx:88
 TGeoEltu.cxx:89
 TGeoEltu.cxx:90
 TGeoEltu.cxx:91
 TGeoEltu.cxx:92
 TGeoEltu.cxx:93
 TGeoEltu.cxx:94
 TGeoEltu.cxx:95
 TGeoEltu.cxx:96
 TGeoEltu.cxx:97
 TGeoEltu.cxx:98
 TGeoEltu.cxx:99
 TGeoEltu.cxx:100
 TGeoEltu.cxx:101
 TGeoEltu.cxx:102
 TGeoEltu.cxx:103
 TGeoEltu.cxx:104
 TGeoEltu.cxx:105
 TGeoEltu.cxx:106
 TGeoEltu.cxx:107
 TGeoEltu.cxx:108
 TGeoEltu.cxx:109
 TGeoEltu.cxx:110
 TGeoEltu.cxx:111
 TGeoEltu.cxx:112
 TGeoEltu.cxx:113
 TGeoEltu.cxx:114
 TGeoEltu.cxx:115
 TGeoEltu.cxx:116
 TGeoEltu.cxx:117
 TGeoEltu.cxx:118
 TGeoEltu.cxx:119
 TGeoEltu.cxx:120
 TGeoEltu.cxx:121
 TGeoEltu.cxx:122
 TGeoEltu.cxx:123
 TGeoEltu.cxx:124
 TGeoEltu.cxx:125
 TGeoEltu.cxx:126
 TGeoEltu.cxx:127
 TGeoEltu.cxx:128
 TGeoEltu.cxx:129
 TGeoEltu.cxx:130
 TGeoEltu.cxx:131
 TGeoEltu.cxx:132
 TGeoEltu.cxx:133
 TGeoEltu.cxx:134
 TGeoEltu.cxx:135
 TGeoEltu.cxx:136
 TGeoEltu.cxx:137
 TGeoEltu.cxx:138
 TGeoEltu.cxx:139
 TGeoEltu.cxx:140
 TGeoEltu.cxx:141
 TGeoEltu.cxx:142
 TGeoEltu.cxx:143
 TGeoEltu.cxx:144
 TGeoEltu.cxx:145
 TGeoEltu.cxx:146
 TGeoEltu.cxx:147
 TGeoEltu.cxx:148
 TGeoEltu.cxx:149
 TGeoEltu.cxx:150
 TGeoEltu.cxx:151
 TGeoEltu.cxx:152
 TGeoEltu.cxx:153
 TGeoEltu.cxx:154
 TGeoEltu.cxx:155
 TGeoEltu.cxx:156
 TGeoEltu.cxx:157
 TGeoEltu.cxx:158
 TGeoEltu.cxx:159
 TGeoEltu.cxx:160
 TGeoEltu.cxx:161
 TGeoEltu.cxx:162
 TGeoEltu.cxx:163
 TGeoEltu.cxx:164
 TGeoEltu.cxx:165
 TGeoEltu.cxx:166
 TGeoEltu.cxx:167
 TGeoEltu.cxx:168
 TGeoEltu.cxx:169
 TGeoEltu.cxx:170
 TGeoEltu.cxx:171
 TGeoEltu.cxx:172
 TGeoEltu.cxx:173
 TGeoEltu.cxx:174
 TGeoEltu.cxx:175
 TGeoEltu.cxx:176
 TGeoEltu.cxx:177
 TGeoEltu.cxx:178
 TGeoEltu.cxx:179
 TGeoEltu.cxx:180
 TGeoEltu.cxx:181
 TGeoEltu.cxx:182
 TGeoEltu.cxx:183
 TGeoEltu.cxx:184
 TGeoEltu.cxx:185
 TGeoEltu.cxx:186
 TGeoEltu.cxx:187
 TGeoEltu.cxx:188
 TGeoEltu.cxx:189
 TGeoEltu.cxx:190
 TGeoEltu.cxx:191
 TGeoEltu.cxx:192
 TGeoEltu.cxx:193
 TGeoEltu.cxx:194
 TGeoEltu.cxx:195
 TGeoEltu.cxx:196
 TGeoEltu.cxx:197
 TGeoEltu.cxx:198
 TGeoEltu.cxx:199
 TGeoEltu.cxx:200
 TGeoEltu.cxx:201
 TGeoEltu.cxx:202
 TGeoEltu.cxx:203
 TGeoEltu.cxx:204
 TGeoEltu.cxx:205
 TGeoEltu.cxx:206
 TGeoEltu.cxx:207
 TGeoEltu.cxx:208
 TGeoEltu.cxx:209
 TGeoEltu.cxx:210
 TGeoEltu.cxx:211
 TGeoEltu.cxx:212
 TGeoEltu.cxx:213
 TGeoEltu.cxx:214
 TGeoEltu.cxx:215
 TGeoEltu.cxx:216
 TGeoEltu.cxx:217
 TGeoEltu.cxx:218
 TGeoEltu.cxx:219
 TGeoEltu.cxx:220
 TGeoEltu.cxx:221
 TGeoEltu.cxx:222
 TGeoEltu.cxx:223
 TGeoEltu.cxx:224
 TGeoEltu.cxx:225
 TGeoEltu.cxx:226
 TGeoEltu.cxx:227
 TGeoEltu.cxx:228
 TGeoEltu.cxx:229
 TGeoEltu.cxx:230
 TGeoEltu.cxx:231
 TGeoEltu.cxx:232
 TGeoEltu.cxx:233
 TGeoEltu.cxx:234
 TGeoEltu.cxx:235
 TGeoEltu.cxx:236
 TGeoEltu.cxx:237
 TGeoEltu.cxx:238
 TGeoEltu.cxx:239
 TGeoEltu.cxx:240
 TGeoEltu.cxx:241
 TGeoEltu.cxx:242
 TGeoEltu.cxx:243
 TGeoEltu.cxx:244
 TGeoEltu.cxx:245
 TGeoEltu.cxx:246
 TGeoEltu.cxx:247
 TGeoEltu.cxx:248
 TGeoEltu.cxx:249
 TGeoEltu.cxx:250
 TGeoEltu.cxx:251
 TGeoEltu.cxx:252
 TGeoEltu.cxx:253
 TGeoEltu.cxx:254
 TGeoEltu.cxx:255
 TGeoEltu.cxx:256
 TGeoEltu.cxx:257
 TGeoEltu.cxx:258
 TGeoEltu.cxx:259
 TGeoEltu.cxx:260
 TGeoEltu.cxx:261
 TGeoEltu.cxx:262
 TGeoEltu.cxx:263
 TGeoEltu.cxx:264
 TGeoEltu.cxx:265
 TGeoEltu.cxx:266
 TGeoEltu.cxx:267
 TGeoEltu.cxx:268
 TGeoEltu.cxx:269
 TGeoEltu.cxx:270
 TGeoEltu.cxx:271
 TGeoEltu.cxx:272
 TGeoEltu.cxx:273
 TGeoEltu.cxx:274
 TGeoEltu.cxx:275
 TGeoEltu.cxx:276
 TGeoEltu.cxx:277
 TGeoEltu.cxx:278
 TGeoEltu.cxx:279
 TGeoEltu.cxx:280
 TGeoEltu.cxx:281
 TGeoEltu.cxx:282
 TGeoEltu.cxx:283
 TGeoEltu.cxx:284
 TGeoEltu.cxx:285
 TGeoEltu.cxx:286
 TGeoEltu.cxx:287
 TGeoEltu.cxx:288
 TGeoEltu.cxx:289
 TGeoEltu.cxx:290
 TGeoEltu.cxx:291
 TGeoEltu.cxx:292
 TGeoEltu.cxx:293
 TGeoEltu.cxx:294
 TGeoEltu.cxx:295
 TGeoEltu.cxx:296
 TGeoEltu.cxx:297
 TGeoEltu.cxx:298
 TGeoEltu.cxx:299
 TGeoEltu.cxx:300
 TGeoEltu.cxx:301
 TGeoEltu.cxx:302
 TGeoEltu.cxx:303
 TGeoEltu.cxx:304
 TGeoEltu.cxx:305
 TGeoEltu.cxx:306
 TGeoEltu.cxx:307
 TGeoEltu.cxx:308
 TGeoEltu.cxx:309
 TGeoEltu.cxx:310
 TGeoEltu.cxx:311
 TGeoEltu.cxx:312
 TGeoEltu.cxx:313
 TGeoEltu.cxx:314
 TGeoEltu.cxx:315
 TGeoEltu.cxx:316
 TGeoEltu.cxx:317
 TGeoEltu.cxx:318
 TGeoEltu.cxx:319
 TGeoEltu.cxx:320
 TGeoEltu.cxx:321
 TGeoEltu.cxx:322
 TGeoEltu.cxx:323
 TGeoEltu.cxx:324
 TGeoEltu.cxx:325
 TGeoEltu.cxx:326
 TGeoEltu.cxx:327
 TGeoEltu.cxx:328
 TGeoEltu.cxx:329
 TGeoEltu.cxx:330
 TGeoEltu.cxx:331
 TGeoEltu.cxx:332
 TGeoEltu.cxx:333
 TGeoEltu.cxx:334
 TGeoEltu.cxx:335
 TGeoEltu.cxx:336
 TGeoEltu.cxx:337
 TGeoEltu.cxx:338
 TGeoEltu.cxx:339
 TGeoEltu.cxx:340
 TGeoEltu.cxx:341
 TGeoEltu.cxx:342
 TGeoEltu.cxx:343
 TGeoEltu.cxx:344
 TGeoEltu.cxx:345
 TGeoEltu.cxx:346
 TGeoEltu.cxx:347
 TGeoEltu.cxx:348
 TGeoEltu.cxx:349
 TGeoEltu.cxx:350
 TGeoEltu.cxx:351
 TGeoEltu.cxx:352
 TGeoEltu.cxx:353
 TGeoEltu.cxx:354
 TGeoEltu.cxx:355
 TGeoEltu.cxx:356
 TGeoEltu.cxx:357
 TGeoEltu.cxx:358
 TGeoEltu.cxx:359
 TGeoEltu.cxx:360
 TGeoEltu.cxx:361
 TGeoEltu.cxx:362
 TGeoEltu.cxx:363
 TGeoEltu.cxx:364
 TGeoEltu.cxx:365
 TGeoEltu.cxx:366
 TGeoEltu.cxx:367
 TGeoEltu.cxx:368
 TGeoEltu.cxx:369
 TGeoEltu.cxx:370
 TGeoEltu.cxx:371
 TGeoEltu.cxx:372
 TGeoEltu.cxx:373
 TGeoEltu.cxx:374
 TGeoEltu.cxx:375
 TGeoEltu.cxx:376
 TGeoEltu.cxx:377
 TGeoEltu.cxx:378
 TGeoEltu.cxx:379
 TGeoEltu.cxx:380
 TGeoEltu.cxx:381
 TGeoEltu.cxx:382
 TGeoEltu.cxx:383
 TGeoEltu.cxx:384
 TGeoEltu.cxx:385
 TGeoEltu.cxx:386
 TGeoEltu.cxx:387
 TGeoEltu.cxx:388
 TGeoEltu.cxx:389
 TGeoEltu.cxx:390
 TGeoEltu.cxx:391
 TGeoEltu.cxx:392
 TGeoEltu.cxx:393
 TGeoEltu.cxx:394
 TGeoEltu.cxx:395
 TGeoEltu.cxx:396
 TGeoEltu.cxx:397
 TGeoEltu.cxx:398
 TGeoEltu.cxx:399
 TGeoEltu.cxx:400
 TGeoEltu.cxx:401
 TGeoEltu.cxx:402
 TGeoEltu.cxx:403
 TGeoEltu.cxx:404
 TGeoEltu.cxx:405
 TGeoEltu.cxx:406
 TGeoEltu.cxx:407
 TGeoEltu.cxx:408
 TGeoEltu.cxx:409
 TGeoEltu.cxx:410
 TGeoEltu.cxx:411
 TGeoEltu.cxx:412
 TGeoEltu.cxx:413
 TGeoEltu.cxx:414
 TGeoEltu.cxx:415
 TGeoEltu.cxx:416
 TGeoEltu.cxx:417
 TGeoEltu.cxx:418
 TGeoEltu.cxx:419
 TGeoEltu.cxx:420
 TGeoEltu.cxx:421
 TGeoEltu.cxx:422
 TGeoEltu.cxx:423
 TGeoEltu.cxx:424
 TGeoEltu.cxx:425
 TGeoEltu.cxx:426
 TGeoEltu.cxx:427
 TGeoEltu.cxx:428
 TGeoEltu.cxx:429
 TGeoEltu.cxx:430
 TGeoEltu.cxx:431
 TGeoEltu.cxx:432
 TGeoEltu.cxx:433
 TGeoEltu.cxx:434
 TGeoEltu.cxx:435
 TGeoEltu.cxx:436
 TGeoEltu.cxx:437
 TGeoEltu.cxx:438
 TGeoEltu.cxx:439
 TGeoEltu.cxx:440
 TGeoEltu.cxx:441
 TGeoEltu.cxx:442
 TGeoEltu.cxx:443
 TGeoEltu.cxx:444
 TGeoEltu.cxx:445
 TGeoEltu.cxx:446
 TGeoEltu.cxx:447
 TGeoEltu.cxx:448
 TGeoEltu.cxx:449
 TGeoEltu.cxx:450
 TGeoEltu.cxx:451
 TGeoEltu.cxx:452
 TGeoEltu.cxx:453
 TGeoEltu.cxx:454
 TGeoEltu.cxx:455
 TGeoEltu.cxx:456
 TGeoEltu.cxx:457
 TGeoEltu.cxx:458
 TGeoEltu.cxx:459
 TGeoEltu.cxx:460
 TGeoEltu.cxx:461
 TGeoEltu.cxx:462
 TGeoEltu.cxx:463
 TGeoEltu.cxx:464
 TGeoEltu.cxx:465
 TGeoEltu.cxx:466
 TGeoEltu.cxx:467
 TGeoEltu.cxx:468
 TGeoEltu.cxx:469
 TGeoEltu.cxx:470
 TGeoEltu.cxx:471
 TGeoEltu.cxx:472
 TGeoEltu.cxx:473
 TGeoEltu.cxx:474
 TGeoEltu.cxx:475
 TGeoEltu.cxx:476
 TGeoEltu.cxx:477
 TGeoEltu.cxx:478
 TGeoEltu.cxx:479
 TGeoEltu.cxx:480
 TGeoEltu.cxx:481
 TGeoEltu.cxx:482
 TGeoEltu.cxx:483
 TGeoEltu.cxx:484
 TGeoEltu.cxx:485
 TGeoEltu.cxx:486
 TGeoEltu.cxx:487
 TGeoEltu.cxx:488
 TGeoEltu.cxx:489
 TGeoEltu.cxx:490
 TGeoEltu.cxx:491
 TGeoEltu.cxx:492
 TGeoEltu.cxx:493
 TGeoEltu.cxx:494
 TGeoEltu.cxx:495
 TGeoEltu.cxx:496
 TGeoEltu.cxx:497
 TGeoEltu.cxx:498
 TGeoEltu.cxx:499
 TGeoEltu.cxx:500
 TGeoEltu.cxx:501
 TGeoEltu.cxx:502
 TGeoEltu.cxx:503
 TGeoEltu.cxx:504
 TGeoEltu.cxx:505
 TGeoEltu.cxx:506
 TGeoEltu.cxx:507
 TGeoEltu.cxx:508
 TGeoEltu.cxx:509
 TGeoEltu.cxx:510
 TGeoEltu.cxx:511
 TGeoEltu.cxx:512
 TGeoEltu.cxx:513
 TGeoEltu.cxx:514
 TGeoEltu.cxx:515
 TGeoEltu.cxx:516
 TGeoEltu.cxx:517
 TGeoEltu.cxx:518
 TGeoEltu.cxx:519
 TGeoEltu.cxx:520
 TGeoEltu.cxx:521
 TGeoEltu.cxx:522
 TGeoEltu.cxx:523
 TGeoEltu.cxx:524
 TGeoEltu.cxx:525
 TGeoEltu.cxx:526
 TGeoEltu.cxx:527
 TGeoEltu.cxx:528
 TGeoEltu.cxx:529
 TGeoEltu.cxx:530
 TGeoEltu.cxx:531
 TGeoEltu.cxx:532
 TGeoEltu.cxx:533
 TGeoEltu.cxx:534
 TGeoEltu.cxx:535
 TGeoEltu.cxx:536
 TGeoEltu.cxx:537
 TGeoEltu.cxx:538
 TGeoEltu.cxx:539
 TGeoEltu.cxx:540
 TGeoEltu.cxx:541
 TGeoEltu.cxx:542
 TGeoEltu.cxx:543
 TGeoEltu.cxx:544
 TGeoEltu.cxx:545
 TGeoEltu.cxx:546
 TGeoEltu.cxx:547
 TGeoEltu.cxx:548
 TGeoEltu.cxx:549
 TGeoEltu.cxx:550
 TGeoEltu.cxx:551
 TGeoEltu.cxx:552
 TGeoEltu.cxx:553
 TGeoEltu.cxx:554