ROOT logo
// @(#)root/geom:$Id: TGeoArb8.cxx 27731 2009-03-09 17:40:56Z brun $
// Author: Andrei Gheata   31/01/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.             *
 *************************************************************************/

#include "Riostream.h"

#include "TGeoManager.h"
#include "TGeoVolume.h"
#include "TGeoArb8.h"
#include "TGeoMatrix.h"
#include "TMath.h"

ClassImp(TGeoArb8)
    
//________________________________________________________________________
// TGeoArb8 - a arbitrary trapezoid with less than 8 vertices standing on
//   two paralel planes perpendicular to Z axis. Parameters :
//            - dz - half length in Z;
//            - xy[8][2] - vector of (x,y) coordinates of vertices
//               - first four points (xy[i][j], i<4, j<2) are the (x,y)
//                 coordinates of the vertices sitting on the -dz plane;
//               - last four points (xy[i][j], i>=4, j<2) are the (x,y)
//                 coordinates of the vertices sitting on the +dz plane;
//   The order of defining the vertices of an arb8 is the following :
//      - point 0 is connected with points 1,3,4
//      - point 1 is connected with points 0,2,5
//      - point 2 is connected with points 1,3,6
//      - point 3 is connected with points 0,2,7
//      - point 4 is connected with points 0,5,7
//      - point 5 is connected with points 1,4,6
//      - point 6 is connected with points 2,5,7
//      - point 7 is connected with points 3,4,6
//   Points can be identical in order to create shapes with less than 
//   8 vertices.
//

//Begin_Html
/*
<img src="gif/t_arb8.gif">
*/
//End_Html

////////////////////////////////////////////////////////////////////////////
//                                                                        //
// TGeoTrap                                                               //
//                                                                        //
// TRAP is a general trapezoid, i.e. one for which the faces perpendicular//
// to z are trapezia and their centres are not the same x, y. It has 11   //
// parameters: the half length in z, the polar angles from the centre of  //
// the face at low z to that at high z, H1 the half length in y at low z, //
// LB1 the half length in x at low z and y low edge, LB2 the half length  //
// in x at low z and y high edge, TH1 the angle w.r.t. the y axis from the//
// centre of low y edge to the centre of the high y edge, and H2, LB2,    //
// LH2, TH2, the corresponding quantities at high z.                      //
//                                                                        //
////////////////////////////////////////////////////////////////////////////
//Begin_Html
/*
<img src="gif/t_trap.gif">
*/
//End_Html
//
//Begin_Html
/*
<img src="gif/t_trapdivZ.gif">
*/
//End_Html

////////////////////////////////////////////////////////////////////////////
//                                                                        //
// TGeoGtra                                                               //
//                                                                        //
// Gtra is a twisted trapezoid, i.e. one for which the faces perpendicular//
// to z are trapezia and their centres are not the same x, y. It has 12   //
// parameters: the half length in z, the polar angles from the centre of  //
// the face at low z to that at high z, twist, H1 the half length in y at low z, //
// LB1 the half length in x at low z and y low edge, LB2 the half length  //
// in x at low z and y high edge, TH1 the angle w.r.t. the y axis from the//
// centre of low y edge to the centre of the high y edge, and H2, LB2,    //
// LH2, TH2, the corresponding quantities at high z.                      //
//                                                                        //
////////////////////////////////////////////////////////////////////////////
//Begin_Html
/*
<img src="gif/t_gtra.gif">
*/
//End_Html
//
//Begin_Html
/*
<img src="gif/t_gtradivstepZ.gif">
*/
//End_Html


//_____________________________________________________________________________
TGeoArb8::TGeoArb8()
{
// Default ctor.
   fDz = 0;
   fTwist = 0;
   for (Int_t i=0; i<8; i++) {
      fXY[i][0] = 0.0;
      fXY[i][1] = 0.0;
   }  
   SetShapeBit(kGeoArb8); 
}

//_____________________________________________________________________________
TGeoArb8::TGeoArb8(Double_t dz, Double_t *vertices)
         :TGeoBBox(0,0,0)
{
// Constructor. If the array of vertices is not null, this should be
// in the format : (x0, y0, x1, y1, ... , x7, y7) 
   fDz = dz;
   fTwist = 0;
   SetShapeBit(kGeoArb8); 
   if (vertices) {
      for (Int_t i=0; i<8; i++) {
         fXY[i][0] = vertices[2*i];
         fXY[i][1] = vertices[2*i+1];
      }
      ComputeTwist();
      ComputeBBox();
   } else {
      for (Int_t i=0; i<8; i++) {
         fXY[i][0] = 0.0;
         fXY[i][1] = 0.0;
      }   
   }
}

//_____________________________________________________________________________
TGeoArb8::TGeoArb8(const char *name, Double_t dz, Double_t *vertices)
         :TGeoBBox(name, 0,0,0)
{
// Named constructor. If the array of vertices is not null, this should be
// in the format : (x0, y0, x1, y1, ... , x7, y7) 
   fDz = dz;
   fTwist = 0;
   SetShapeBit(kGeoArb8); 
   if (vertices) {
      for (Int_t i=0; i<8; i++) {
         fXY[i][0] = vertices[2*i];
         fXY[i][1] = vertices[2*i+1];
      }
      ComputeTwist();
      ComputeBBox();
   } else {
      for (Int_t i=0; i<8; i++) {
         fXY[i][0] = 0.0;
         fXY[i][1] = 0.0;
      }   
   }
}

//_____________________________________________________________________________
TGeoArb8::TGeoArb8(const TGeoArb8& ga8) :
  TGeoBBox(ga8),
  fDz(ga8.fDz),
  fTwist(ga8.fTwist)
{
   //copy constructor
   for(Int_t i=0; i<8; i++) {
      fXY[i][0]=ga8.fXY[i][0];
      fXY[i][1]=ga8.fXY[i][1];
   }
}

//_____________________________________________________________________________
TGeoArb8& TGeoArb8::operator=(const TGeoArb8& ga8) 
{
   //assignment operator
   if(this!=&ga8) {
      TGeoBBox::operator=(ga8);
      fDz=ga8.fDz;
      fTwist=ga8.fTwist;
      for(Int_t i=0; i<8; i++) {
         fXY[i][0]=ga8.fXY[i][0];
         fXY[i][1]=ga8.fXY[i][1];
      }
   } 
   return *this;
}

//_____________________________________________________________________________
TGeoArb8::~TGeoArb8()
{
// Destructor.
   if (fTwist) delete [] fTwist;
}

//_____________________________________________________________________________
Double_t TGeoArb8::Capacity() const
{
// Computes capacity of the shape in [length^3].
   Int_t i,j;
   Double_t capacity = 0;
   for (i=0; i<4; i++) {
      j = (i+1)%4;
      capacity += 0.25*fDz*((fXY[i][0]+fXY[i+4][0])*(fXY[j][1]+fXY[j+4][1]) -
                            (fXY[j][0]+fXY[j+4][0])*(fXY[i][1]+fXY[i+4][1]) +
                    (1./3)*((fXY[i+4][0]-fXY[i][0])*(fXY[j+4][1]-fXY[j][1]) -
                            (fXY[j][0]-fXY[j+4][0])*(fXY[i][1]-fXY[i+4][1])));
   }
   return TMath::Abs(capacity);
}                                

//_____________________________________________________________________________
void TGeoArb8::ComputeBBox()
{
// Computes bounding box for an Arb8 shape.
   Double_t xmin, xmax, ymin, ymax;
   xmin = xmax = fXY[0][0];
   ymin = ymax = fXY[0][1];
   
   for (Int_t i=1; i<8; i++) {
      if (xmin>fXY[i][0]) xmin=fXY[i][0];
      if (xmax<fXY[i][0]) xmax=fXY[i][0];
      if (ymin>fXY[i][1]) ymin=fXY[i][1];
      if (ymax<fXY[i][1]) ymax=fXY[i][1];
   }
   fDX = 0.5*(xmax-xmin);
   fDY = 0.5*(ymax-ymin);
   fDZ = fDz;
   fOrigin[0] = 0.5*(xmax+xmin);
   fOrigin[1] = 0.5*(ymax+ymin);
   fOrigin[2] = 0;
   SetShapeBit(kGeoClosedShape);
}   

//_____________________________________________________________________________
void TGeoArb8::ComputeTwist()
{
// Computes tangents of twist angles (angles between projections on XY plane
// of corresponding -dz +dz edges). Computes also if the vertices are defined
// clockwise or anti-clockwise.
   Double_t twist[4];
   Bool_t twisted = kFALSE;
   Double_t dx1, dy1, dx2, dy2;
   for (Int_t i=0; i<4; i++) {
      dx1 = fXY[(i+1)%4][0]-fXY[i][0];
      dy1 = fXY[(i+1)%4][1]-fXY[i][1];
      if (TMath::Abs(dx1)<TGeoShape::Tolerance() && TMath::Abs(dy1)<TGeoShape::Tolerance()) {
         twist[i] = 0;
         continue;
      }   
      dx2 = fXY[4+(i+1)%4][0]-fXY[4+i][0];
      dy2 = fXY[4+(i+1)%4][1]-fXY[4+i][1];
      if (TMath::Abs(dx2)<TGeoShape::Tolerance() && TMath::Abs(dy2)<TGeoShape::Tolerance()) {
         twist[i] = 0;
         continue;
      }
      twist[i] = dy1*dx2 - dx1*dy2;
      if (TMath::Abs(twist[i])<TGeoShape::Tolerance()) {
         twist[i] = 0;
         continue;
      }
      twist[i] = TMath::Sign(1.,twist[i]);
      twisted = kTRUE;
   }
   if (twisted) {
      if (fTwist) delete [] fTwist;
      fTwist = new Double_t[4];
      memcpy(fTwist, &twist[0], 4*sizeof(Double_t));
   }   
   Double_t sum1 = 0.;
   Double_t sum2 = 0.;
   Int_t j;
   for (Int_t i=0; i<4; i++) {
      j = (i+1)%4;
      sum1 += fXY[i][0]*fXY[j][1]-fXY[j][0]*fXY[i][1];
      sum2 += fXY[i+4][0]*fXY[j+4][1]-fXY[j+4][0]*fXY[i+4][1];
   }
   if (sum1*sum2 < -TGeoShape::Tolerance()) {
      Fatal("ComputeTwist", "Shape %s type Arb8: Lower/upper faces defined with opposite clockwise", GetName());
      return;
   }
   if (sum1>0.) {
      Error("ComputeTwist", "Shape %s type Arb8: Vertices must be defined clockwise in XY planes. Re-ordering...", GetName());
      Double_t xtemp, ytemp;
      xtemp = fXY[1][0];
      ytemp = fXY[1][1];
      fXY[1][0] = fXY[3][0];
      fXY[1][1] = fXY[3][1];
      fXY[3][0] = xtemp;
      fXY[3][1] = ytemp;
      xtemp = fXY[5][0];
      ytemp = fXY[5][1];
      fXY[5][0] = fXY[7][0];
      fXY[5][1] = fXY[7][1];
      fXY[7][0] = xtemp;
      fXY[7][1] = ytemp;
   } 
   // Check for illegal crossings.
   Bool_t illegal_cross = kFALSE;
   illegal_cross = TGeoShape::IsSegCrossing(fXY[0][0],fXY[0][1],fXY[1][0],fXY[1][1],
                                            fXY[2][0],fXY[2][1],fXY[3][0],fXY[3][1]);
   if (!illegal_cross) 
   illegal_cross = TGeoShape::IsSegCrossing(fXY[4][0],fXY[4][1],fXY[5][0],fXY[5][1],
                                            fXY[6][0],fXY[6][1],fXY[7][0],fXY[7][1]);
   if (illegal_cross) {
      Error("ComputeTwist", "Shape %s type Arb8: Malformed polygon with crossing opposite segments", GetName());
      InspectShape();
   }
}

//_____________________________________________________________________________
Double_t TGeoArb8::GetTwist(Int_t iseg) const
{
// Get twist for segment I in range [0,3]
   if (!fTwist) return 0.;
   if (iseg<0 || iseg>3) return 0.;
   return fTwist[iseg];
}   

//_____________________________________________________________________________
void TGeoArb8::ComputeNormal(Double_t *point, Double_t *dir, Double_t *norm)
{
// Compute normal to closest surface from POINT. 
   Double_t safe = TGeoShape::Big();
   Double_t safc;
   Int_t i;          // current facette index
   Double_t x0, y0, z0, x1, y1, z1, x2, y2;
   Double_t ax, ay, az, bx, by;
   Double_t fn;
   safc = fDz-TMath::Abs(point[2]);
   if (safc<1E-4) {
      memset(norm,0,3*sizeof(Double_t));
      norm[2] = (dir[2]>0)?1:(-1);
      return;
   }   
   Double_t vert[8], lnorm[3];
   SetPlaneVertices(point[2], vert);
   //---> compute safety for lateral planes
   for (i=0; i<4; i++) {
      x0 = vert[2*i];
      y0 = vert[2*i+1];
      z0 = point[2];
      x1 = fXY[i+4][0];
      y1 = fXY[i+4][1];
      z1 = fDz;
      ax = x1-x0;
      ay = y1-y0;
      az = z1-z0;
      x2 = vert[2*((i+1)%4)];
      y2 = vert[2*((i+1)%4)+1];
      bx = x2-x0;
      by = y2-y0;

      lnorm[0] = -az*by;
      lnorm[1] = az*bx;
      lnorm[2] = ax*by-ay*bx;
      fn = TMath::Sqrt(lnorm[0]*lnorm[0]+lnorm[1]*lnorm[1]+lnorm[2]*lnorm[2]);
      if (fn<1E-10) continue;
      lnorm[0] /= fn;
      lnorm[1] /= fn;
      lnorm[2] /= fn;
      safc = (x0-point[0])*lnorm[0]+(y0-point[1])*lnorm[1]+(z0-point[2])*lnorm[2];
      safc = TMath::Abs(safc);
//      printf("plane %i : (%g, %g, %g) safe=%g\n", i, lnorm[0],lnorm[1],lnorm[2],safc);
      if (safc<safe) {
         safe = safc;
         memcpy(norm,lnorm,3*sizeof(Double_t));
      }      
   }
   if (dir[0]*norm[0]+dir[1]*norm[1]+dir[2]*norm[2] < 0) {
      norm[0] = -norm[0];
      norm[1] = -norm[1];
      norm[2] = -norm[2];
   }
}

//_____________________________________________________________________________
Bool_t TGeoArb8::Contains(Double_t *point) const
{
// Test if point is inside this shape.
   // first check Z range
   if (TMath::Abs(point[2]) > fDz) return kFALSE;
   // compute intersection between Z plane containing point and the arb8
   Double_t poly[8];
//   memset(&poly[0], 0, 8*sizeof(Double_t));
   //SetPlaneVertices(point[2], &poly[0]);
   Double_t cf = 0.5*(fDz-point[2])/fDz;
   Int_t i;
   for (i=0; i<4; i++) {
      poly[2*i]   = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
      poly[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
   }
   return InsidePolygon(point[0],point[1],poly);
}

//_____________________________________________________________________________
Double_t TGeoArb8::DistToPlane(Double_t *point, Double_t *dir, Int_t ipl, Bool_t in) const 
{
// Computes distance to plane ipl :
// ipl=0 : points 0,4,1,5
// ipl=1 : points 1,5,2,6
// ipl=2 : points 2,6,3,7
// ipl=3 : points 3,7,0,4
   Double_t xa,xb,xc,xd;
   Double_t ya,yb,yc,yd;
   Double_t eps = 100.*TGeoShape::Tolerance();
   Int_t j = (ipl+1)%4;
   xa=fXY[ipl][0];
   ya=fXY[ipl][1];
   xb=fXY[ipl+4][0];
   yb=fXY[ipl+4][1];
   xc=fXY[j][0];
   yc=fXY[j][1];
   xd=fXY[4+j][0];
   yd=fXY[4+j][1];
   Double_t dz2 =0.5/fDz;
   Double_t tx1 =dz2*(xb-xa);
   Double_t ty1 =dz2*(yb-ya);
   Double_t tx2 =dz2*(xd-xc);
   Double_t ty2 =dz2*(yd-yc);
   Double_t dzp =fDz+point[2];
   Double_t xs1 =xa+tx1*dzp;
   Double_t ys1 =ya+ty1*dzp;
   Double_t xs2 =xc+tx2*dzp;
   Double_t ys2 =yc+ty2*dzp;
   Double_t dxs =xs2-xs1;
   Double_t dys =ys2-ys1;
   Double_t dtx =tx2-tx1;
   Double_t dty =ty2-ty1;
   Double_t a=(dtx*dir[1]-dty*dir[0]+(tx1*ty2-tx2*ty1)*dir[2])*dir[2];
   Double_t b=dxs*dir[1]-dys*dir[0]+(dtx*point[1]-dty*point[0]+ty2*xs1-ty1*xs2
              +tx1*ys2-tx2*ys1)*dir[2];
   Double_t c=dxs*point[1]-dys*point[0]+xs1*ys2-xs2*ys1;
   Double_t s=TGeoShape::Big();
   Double_t x1,x2,y1,y2,xp,yp,zi;
   if (TMath::Abs(a)<eps) {           
      if (TMath::Abs(b)<eps) return TGeoShape::Big();
      s=-c/b;
      if (s>eps) {
         if (in) return s;
         zi=point[2]+s*dir[2];
         if (TMath::Abs(zi)<fDz) {
            x1=xs1+tx1*dir[2]*s;
            x2=xs2+tx2*dir[2]*s;
            xp=point[0]+s*dir[0];
            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 TGeoShape::Big();
   }      
   Double_t d=b*b-4*a*c;
   if (d>=0) {
      if (a>0) s=0.5*(-b-TMath::Sqrt(d))/a;
      else     s=0.5*(-b+TMath::Sqrt(d))/a;
      if (s>eps) {
         if (in) return s;
         zi=point[2]+s*dir[2];
         if (TMath::Abs(zi)<fDz) {
            x1=xs1+tx1*dir[2]*s;
            x2=xs2+tx2*dir[2]*s;
            xp=point[0]+s*dir[0];
            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;
         }      
      }
      if (a>0) s=0.5*(-b+TMath::Sqrt(d))/a;
      else     s=0.5*(-b-TMath::Sqrt(d))/a;
      if (s>eps) {
         if (in) return s;
         zi=point[2]+s*dir[2];
         if (TMath::Abs(zi)<fDz) {
            x1=xs1+tx1*dir[2]*s;
            x2=xs2+tx2*dir[2]*s;
            xp=point[0]+s*dir[0];
            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 TGeoShape::Big();
}      
      
//_____________________________________________________________________________
Double_t TGeoArb8::DistFromOutside(Double_t *point, Double_t *dir, Int_t /*iact*/, Double_t step, Double_t * /*safe*/) const
{
// Computes distance from outside point to surface of the shape.
   Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
   if (sdist>=step) return TGeoShape::Big();
   Double_t dist[5];
   // check lateral faces
   Int_t i;
   for (i=0; i<4; i++) {
      dist[i]=DistToPlane(point, dir, i, kFALSE);  
   }   
   // check Z planes
   dist[4]=TGeoShape::Big();
   if (TMath::Abs(point[2])>fDz) {
      if (!TGeoShape::IsSameWithinTolerance(dir[2],0)) {
         Double_t pt[3];
         if (point[2]>0) {
            dist[4] = (fDz-point[2])/dir[2];
            pt[2]=fDz;
         } else {   
            dist[4] = (-fDz-point[2])/dir[2];
            pt[2]=-fDz;
         }   
         if (dist[4]<0) {
            dist[4]=TGeoShape::Big();
         } else {   
            for (Int_t j=0; j<2; j++) pt[j]=point[j]+dist[4]*dir[j];
            if (!Contains(&pt[0])) dist[4]=TGeoShape::Big();
         }   
      }
   }   
   Double_t distmin = dist[0];
   for (i=1;i<5;i++) if (dist[i] < distmin) distmin = dist[i];
   return distmin;
}   

//_____________________________________________________________________________
Double_t TGeoArb8::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 shape.
#ifdef OLDALGORITHM
   Int_t i;
   Double_t dist[6];
   dist[0]=dist[1]=TGeoShape::Big();
   if (dir[2]<0) {
      dist[0]=(-fDz-point[2])/dir[2];
   } else {
      if (dir[2]>0) dist[1]=(fDz-point[2])/dir[2];
   }      
   for (i=0; i<4; i++) {
      dist[i+2]=DistToPlane(point, dir, i, kTRUE);   
   }   
      
   Double_t distmin = dist[0];
   for (i=1;i<6;i++) if (dist[i] < distmin) distmin = dist[i];
   return distmin;
#else
// compute distance to plane ipl :
// ipl=0 : points 0,4,1,5
// ipl=1 : points 1,5,2,6
// ipl=2 : points 2,6,3,7
// ipl=3 : points 3,7,0,4
   Double_t distmin;
   Bool_t lateral_cross = kFALSE;
   if (dir[2]<0) {
      distmin=(-fDz-point[2])/dir[2];
   } else {
      if (dir[2]>0) distmin =(fDz-point[2])/dir[2];
      else          distmin = TGeoShape::Big();
   }      
   Double_t dz2 =0.5/fDz;
   Double_t xa,xb,xc,xd;
   Double_t ya,yb,yc,yd;
   Double_t eps = 100.*TGeoShape::Tolerance();
   for (Int_t ipl=0;ipl<4;ipl++) {
      Int_t j = (ipl+1)%4;
      xa=fXY[ipl][0];
      ya=fXY[ipl][1];
      xb=fXY[ipl+4][0];
      yb=fXY[ipl+4][1];
      xc=fXY[j][0];
      yc=fXY[j][1];
      xd=fXY[4+j][0];
      yd=fXY[4+j][1];
      
      Double_t tx1 =dz2*(xb-xa);
      Double_t ty1 =dz2*(yb-ya);
      Double_t tx2 =dz2*(xd-xc);
      Double_t ty2 =dz2*(yd-yc);
      Double_t dzp =fDz+point[2];
      Double_t xs1 =xa+tx1*dzp;
      Double_t ys1 =ya+ty1*dzp;
      Double_t xs2 =xc+tx2*dzp;
      Double_t ys2 =yc+ty2*dzp;
      Double_t dxs =xs2-xs1;
      Double_t dys =ys2-ys1;
      Double_t dtx =tx2-tx1;
      Double_t dty =ty2-ty1;
      Double_t a=(dtx*dir[1]-dty*dir[0]+(tx1*ty2-tx2*ty1)*dir[2])*dir[2];
      Double_t b=dxs*dir[1]-dys*dir[0]+(dtx*point[1]-dty*point[0]+ty2*xs1-ty1*xs2
              +tx1*ys2-tx2*ys1)*dir[2];
      Double_t c=dxs*point[1]-dys*point[0]+xs1*ys2-xs2*ys1;
      Double_t s=TGeoShape::Big();
      if (TMath::Abs(a)<eps) {           
         if (TMath::Abs(b)<eps) continue;
         s=-c/b;
         if (s>eps && s < distmin) {
            distmin =s;
            lateral_cross=kTRUE;
         }   
         continue;
      }
      Double_t d=b*b-4*a*c;
      if (d>=0.) {
         if (a>0) s=0.5*(-b-TMath::Sqrt(d))/a;
         else     s=0.5*(-b+TMath::Sqrt(d))/a;
         if (s>eps) {
            if (s < distmin) {
               distmin = s;
               lateral_cross = kTRUE;
            }   
         } else {
            if (a>0) s=0.5*(-b+TMath::Sqrt(d))/a;
            else     s=0.5*(-b-TMath::Sqrt(d))/a;
            if (s>eps && s < distmin) {
               distmin =s;
               lateral_cross = kTRUE;
            }   
         }
      }
   }

   if (!lateral_cross) {
      // We have to make sure that track crosses the top or bottom.
      if (distmin > 1.E10) return TGeoShape::Tolerance();
      Double_t pt[2];
      pt[0] = point[0]+distmin*dir[0];
      pt[1] = point[1]+distmin*dir[1];
      // Check if propagated point is in the polygon
      Double_t poly[8];
      Int_t i = 0;
      if (dir[2]>0.) i=4;
      for (Int_t j=0; j<4; j++) {
         poly[2*j]   = fXY[j+i][0];
         poly[2*j+1] = fXY[j+i][1];
      }
      if (!InsidePolygon(pt[0],pt[1],poly)) return TGeoShape::Tolerance();
   }   
   return distmin;
#endif
}   

//_____________________________________________________________________________
TGeoVolume *TGeoArb8::Divide(TGeoVolume *voldiv, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/, 
                             Double_t /*start*/, Double_t /*step*/) 
{
// Divide this shape along one axis.
   Error("Divide", "Division of an arbitrary trapezoid not implemented");
   return voldiv;
}      

//_____________________________________________________________________________
Double_t TGeoArb8::GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
{
// Get shape range on a given axis.
   xlo = 0;
   xhi = 0;
   Double_t dx = 0;
   if (iaxis==3) {
      xlo = -fDz;
      xhi = fDz;
      dx = xhi-xlo;
      return dx;
   }
   return dx;
}
      
//_____________________________________________________________________________
void TGeoArb8::GetBoundingCylinder(Double_t *param) const
{
//--- Fill vector param[4] with the bounding cylinder parameters. The order
// is the following : Rmin, Rmax, Phi1, Phi2
   //--- first compute rmin/rmax
   Double_t rmaxsq = 0;
   Double_t rsq;
   Int_t i;
   for (i=0; i<8; i++) {
      rsq = fXY[i][0]*fXY[i][0] + fXY[i][1]*fXY[i][1];
      rmaxsq = TMath::Max(rsq, rmaxsq);
   }
   param[0] = 0.;                  // Rmin
   param[1] = rmaxsq;              // Rmax
   param[2] = 0.;                  // Phi1
   param[3] = 360.;                // Phi2 
}   

//_____________________________________________________________________________
Int_t TGeoArb8::GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const
{
// Fills real parameters of a positioned box inside this arb8. Returns 0 if successfull.
   dx=dy=dz=0;
   if (mat->IsRotation()) {
      Error("GetFittingBox", "cannot handle parametrized rotated volumes");
      return 1; // ### rotation not accepted ###
   }
   //--> translate the origin of the parametrized box to the frame of this box.
   Double_t origin[3];
   mat->LocalToMaster(parambox->GetOrigin(), origin);
   if (!Contains(origin)) {
      Error("GetFittingBox", "wrong matrix - parametrized box is outside this");
      return 1; // ### wrong matrix ###
   }
   //--> now we have to get the valid range for all parametrized axis
   Double_t dd[3];
   dd[0] = parambox->GetDX();
   dd[1] = parambox->GetDY();
   dd[2] = parambox->GetDZ();
   //-> check if Z range is fixed
   if (dd[2]<0) {
      dd[2] = TMath::Min(origin[2]+fDz, fDz-origin[2]); 
      if (dd[2]<0) {
         Error("GetFittingBox", "wrong matrix");
         return 1;
      }
   }
   if (dd[0]>=0 && dd[1]>=0) {
      dx = dd[0];
      dy = dd[1];
      dz = dd[2];
      return 0;
   }
   //-> check now vertices at Z = origin[2] +/- dd[2]
   Double_t upper[8];
   Double_t lower[8];
   SetPlaneVertices(origin[2]-dd[2], lower);
   SetPlaneVertices(origin[2]+dd[2], upper);
   Double_t ddmin=TGeoShape::Big();
   for (Int_t iaxis=0; iaxis<2; iaxis++) {
      if (dd[iaxis]>=0) continue;
      ddmin=TGeoShape::Big();
      for (Int_t ivert=0; ivert<4; ivert++) {
         ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-lower[2*ivert+iaxis]));
         ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-upper[2*ivert+iaxis]));
      }
      dd[iaxis] = ddmin;
   }
   dx = dd[0];
   dy = dd[1];
   dz = dd[2];
   return 0;
}   

//_____________________________________________________________________________
void TGeoArb8::GetPlaneNormal(Double_t *p1, Double_t *p2, Double_t *p3, Double_t *norm)
{
// Computes normal to plane defined by P1, P2 and P3
   Double_t cross = 0.;
   Double_t v1[3], v2[3];
   Int_t i;
   for (i=0; i<3; i++) {
      v1[i] = p2[i] - p1[i];
      v2[i] = p3[i] - p1[i];
   }
   norm[0] = v1[1]*v2[2]-v1[2]*v2[1];
   cross += norm[0]*norm[0];
   norm[1] = v1[2]*v2[0]-v1[0]*v2[2];
   cross += norm[1]*norm[1];
   norm[2] = v1[0]*v2[1]-v1[1]*v2[0];
   cross += norm[2]*norm[2];
   if (TMath::Abs(cross) < TGeoShape::Tolerance()) return;
   cross = 1./TMath::Sqrt(cross);
   for (i=0; i<3; i++) norm[i] *= cross;
}   

//_____________________________________________________________________________
Bool_t TGeoArb8::GetPointsOnFacet(Int_t /*index*/, Int_t /*npoints*/, Double_t * /* array */) const
{
// Fills array with n random points located on the surface of indexed facet.
// The output array must be provided with a length of minimum 3*npoints. Returns
// true if operation succeeded.
// Possible index values:
//    0 - all facets togeather
//    1 to 6 - facet index from bottom to top Z
   return kFALSE;
/*
   if (index<0 || index>6) return kFALSE;
   if (index==0) {
   // Just generate same number of points on each facet
      Int_t npts = npoints/6.;
      Int_t count = 0;
      for (Int_t ifacet=0; ifacet<6; ifacet++) {
         if (GetPointsOnFacet(ifacet+1, npts, &array[3*count])) count += npts;
         if (ifacet<5) npts = (npoints-count)/(5.-ifacet);
      }   
      if (count>0) return kTRUE;
      return kFALSE;
   }   
   Double_t z, cf;
   Double_t xmin=TGeoShape::Big();
   Double_t xmax=-xmin;
   Double_t ymin=TGeoShape::Big();
   Double_t ymax=-ymin;
   Double_t dy=0.;
   Double_t poly[8];
   Double_t point[2];
   Int_t i;
   if (index==1 || index==6) {
      z = (index==1)?-fDz:fDz;
      cf = 0.5*(fDz-z)/fDz;
      for (i=0; i<4; i++) {
         poly[2*i]   = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
         poly[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
         xmin = TMath::Min(xmin, poly[2*i]);
         xmax = TMath::Max(xmax, poly[2*i]);
         ymin = TMath::Min(ymin, poly[2*i]);
         ymax = TMath::Max(ymax, poly[2*i]);
      }
   }
   Int_t nshoot = 0;
   Int_t nmiss = 0;
   for (i=0; i<npoints; i++) {
      Double_t *point = &array[3*i];
      switch (surfindex) {
         case 1:
         case 6:
            while (nmiss<1000) {
               point[0] = xmin + (xmax-xmin)*gRandom->Rndm();
               point[1] = ymin + (ymax-ymin)*gRandom->Rndm();
            }   

   return InsidePolygon(point[0],point[1],poly);
*/
}   

//_____________________________________________________________________________
Bool_t TGeoArb8::InsidePolygon(Double_t x, Double_t y, Double_t *pts)
{
// Finds if a point in XY plane is inside the polygon defines by PTS.
   Int_t i,j;
   Double_t x1,y1,x2,y2;
   Double_t cross;
   for (i=0; i<4; i++) {
      j = (i+1)%4;
      x1 = pts[i<<1];
      y1 = pts[(i<<1)+1];
      x2 = pts[j<<1];
      y2 = pts[(j<<1)+1];
      cross = (x-x1)*(y2-y1)-(y-y1)*(x2-x1);
      if (cross<0) return kFALSE;
   }
   return kTRUE;   
}

//_____________________________________________________________________________
void TGeoArb8::InspectShape() const
{
// Prints shape parameters
   printf("*** Shape %s: TGeoArb8 ***\n", GetName());
   if (IsTwisted()) printf("  = TWISTED\n");
   for (Int_t ip=0; ip<8; ip++) {
      printf("    point #%i : x=%11.5f y=%11.5f z=%11.5f\n", 
             ip, fXY[ip][0], fXY[ip][1], fDz*((ip<4)?-1:1));
   }
   printf(" Bounding box:\n");
   TGeoBBox::InspectShape();
}

//_____________________________________________________________________________
Double_t TGeoArb8::Safety(Double_t *point, Bool_t in) const
{
// Computes the closest distance from given point to this shape.
   Double_t safz = fDz-TMath::Abs(point[2]);
   if (!in) safz = -safz;
   Int_t iseg;
   Double_t safe = TGeoShape::Big();
   Double_t lsq, ssq, dx, dy, dpx, dpy, u;
   if (IsTwisted()) {
      if (!in) {
         if (!TGeoBBox::Contains(point)) return TGeoBBox::Safety(point,kFALSE);
      }
      // Point is also in the bounding box ;-(  
      // Compute closest distance to any segment
      Double_t vert[8];
      Double_t *p1, *p2;
      Int_t isegmin=0;
      Double_t umin = 0.;
      SetPlaneVertices (point[2], vert);
      for (iseg=0; iseg<4; iseg++) {
         if (safe<TGeoShape::Tolerance()) return 0.;
         p1 = &vert[2*iseg];
         p2 = &vert[2*((iseg+1)%4)];
         dx = p2[0] - p1[0];
         dy = p2[1] - p1[1];
         dpx = point[0] - p1[0];
         dpy = point[1] - p1[1];
      
         lsq = dx*dx + dy*dy;
         u = (dpx*dx + dpy*dy)/lsq;
         if (u>1) {
            dpx = point[0]-p2[0];
            dpy = point[1]-p2[1];
         } else {
            if (u>=0) {
               dpx -= u*dx;
               dpy -= u*dy;
            }
         }
         ssq = dpx*dpx + dpy*dpy;      
         if (ssq < safe) {
            isegmin = iseg;
            umin = u;
            safe = ssq;
         }   
      }
      if (umin<0) umin = 0.;
      if (umin>1) {
         isegmin = (isegmin+1)%4;
         umin = 0.;
      }
      Int_t i1 = isegmin;
      Int_t i2 = (isegmin+1)%4;
      Double_t dx1 = fXY[i2][0]-fXY[i1][0];   
      Double_t dx2 = fXY[i2+4][0]-fXY[i1+4][0];   
      Double_t dy1 = fXY[i2][1]-fXY[i1][1];   
      Double_t dy2 = fXY[i2+4][1]-fXY[i1+4][1];
      dx = dx1 + umin*(dx2-dx1);
      dy = dy1 + umin*(dy2-dy1);
      safe *= 1.- 4.*fDz*fDz/(dx*dx+dy*dy+4.*fDz*fDz);
      safe = TMath::Sqrt(safe);      
      return safe;   
   }  
      
   Double_t saf[5];
   saf[0] = safz;
      
   for (iseg=0; iseg<4; iseg++) saf[iseg+1] = SafetyToFace(point,iseg,in);
   if (in) safe = saf[TMath::LocMin(5, saf)];
   else    safe = saf[TMath::LocMax(5, saf)];   
   if (safe<0) return 0.;
   return safe;
}

//_____________________________________________________________________________
Double_t TGeoArb8::SafetyToFace(Double_t *point, Int_t iseg, Bool_t in) const
{
// Estimate safety to lateral plane defined by segment iseg in range [0,3]
// Might be negative: plane seen only from inside.
   Double_t vertices[12];
   Int_t ipln = (iseg+1)%4;
   // point 1
   vertices[0] = fXY[iseg][0];
   vertices[1] = fXY[iseg][1];
   vertices[2] = -fDz;
   // point 2
   vertices[3] = fXY[ipln][0];
   vertices[4] = fXY[ipln][1];
   vertices[5] = -fDz;
   // point 3
   vertices[6] = fXY[ipln+4][0];
   vertices[7] = fXY[ipln+4][1];
   vertices[8] = fDz;
   // point 4
   vertices[9] = fXY[iseg+4][0];
   vertices[10] = fXY[iseg+4][1];
   vertices[11] = fDz;
   Double_t safe;
   Double_t norm[3];
   Double_t *p1, *p2, *p3;
   p1 = &vertices[0];
   p2 = &vertices[9];
   p3 = &vertices[6];
   if (IsSamePoint(p2,p3)) {
      p3 = &vertices[3];
      if (IsSamePoint(p1,p3)) return -TGeoShape::Big(); // skip single segment
   }
   GetPlaneNormal(p1,p2,p3,norm);
   safe = (point[0]-p1[0])*norm[0]+(point[1]-p1[1])*norm[1]+(point[2]-p1[2])*norm[2];
   if (in) return (-safe);
   return safe;
}
   
//_____________________________________________________________________________
void TGeoArb8::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 << "   dz       = " << fDz << ";" << endl;
   out << "   vert[0]  = " << fXY[0][0] << ";" << endl;
   out << "   vert[1]  = " << fXY[0][1] << ";" << endl;
   out << "   vert[2]  = " << fXY[1][0] << ";" << endl;
   out << "   vert[3]  = " << fXY[1][1] << ";" << endl;
   out << "   vert[4]  = " << fXY[2][0] << ";" << endl;
   out << "   vert[5]  = " << fXY[2][1] << ";" << endl;
   out << "   vert[6]  = " << fXY[3][0] << ";" << endl;
   out << "   vert[7]  = " << fXY[3][1] << ";" << endl;
   out << "   vert[8]  = " << fXY[4][0] << ";" << endl;
   out << "   vert[9]  = " << fXY[4][1] << ";" << endl;
   out << "   vert[10] = " << fXY[5][0] << ";" << endl;
   out << "   vert[11] = " << fXY[5][1] << ";" << endl;
   out << "   vert[12] = " << fXY[6][0] << ";" << endl;
   out << "   vert[13] = " << fXY[6][1] << ";" << endl;
   out << "   vert[14] = " << fXY[7][0] << ";" << endl;
   out << "   vert[15] = " << fXY[7][1] << ";" << endl;
   out << "   TGeoShape *" << GetPointerName() << " = new TGeoArb8(\"" << GetName() << "\", dz,vert);" << endl;
   TObject::SetBit(TGeoShape::kGeoSavePrimitive);
}

//_____________________________________________________________________________
void TGeoArb8::SetPlaneVertices(Double_t zpl, Double_t *vertices) const
{
 // Computes intersection points between plane at zpl and non-horizontal edges.
   Double_t cf = 0.5*(fDz-zpl)/fDz;
   for (Int_t i=0; i<4; i++) {
      vertices[2*i]   = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
      vertices[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
   }
}

//_____________________________________________________________________________
void TGeoArb8::SetDimensions(Double_t *param)
{
// Set all arb8 params in one step.
// param[0] = dz
// param[1] = x0
// param[2] = y0
// ...
   fDz      = param[0];
   for (Int_t i=0; i<8; i++) {
      fXY[i][0] = param[2*i+1];
      fXY[i][1] = param[2*i+2];
   }
   ComputeTwist();
   ComputeBBox();
}   

//_____________________________________________________________________________
void TGeoArb8::SetPoints(Double_t *points) const
{
// Creates arb8 mesh points
   for (Int_t i=0; i<8; i++) {
      points[3*i] = fXY[i][0];
      points[3*i+1] = fXY[i][1];
      points[3*i+2] = (i<4)?-fDz:fDz;
   }
}

//_____________________________________________________________________________
void TGeoArb8::SetPoints(Float_t *points) const
{
// Creates arb8 mesh points
   for (Int_t i=0; i<8; i++) {
      points[3*i] = fXY[i][0];
      points[3*i+1] = fXY[i][1];
      points[3*i+2] = (i<4)?-fDz:fDz;
   }
}

//_____________________________________________________________________________
void TGeoArb8::SetVertex(Int_t vnum, Double_t x, Double_t y)
{
//  Set values for a given vertex.
   if (vnum<0 || vnum >7) {
      Error("SetVertex", "Invalid vertex number");
      return;
   }
   fXY[vnum][0] = x;
   fXY[vnum][1] = y;
   if (vnum == 7) {
      ComputeTwist();
      ComputeBBox();
   }
}

//_____________________________________________________________________________
void TGeoArb8::Sizeof3D() const
{
// Fill size of this 3-D object
   TGeoBBox::Sizeof3D();
}

//_____________________________________________________________________________
void TGeoArb8::Streamer(TBuffer &R__b)
{
   // Stream an object of class TGeoManager.
   if (R__b.IsReading()) {
      R__b.ReadClassBuffer(TGeoArb8::Class(), this);
      ComputeTwist();
   } else {
      R__b.WriteClassBuffer(TGeoArb8::Class(), this);
   }
}   

ClassImp(TGeoTrap)

//_____________________________________________________________________________
TGeoTrap::TGeoTrap()
{
// Default ctor
   fDz = 0;
   fTheta = 0;
   fPhi = 0;
   fH1 = fH2 = fBl1 = fBl2 = fTl1 = fTl2 = fAlpha1 = fAlpha2 = 0;
}

//_____________________________________________________________________________
TGeoTrap::TGeoTrap(Double_t dz, Double_t theta, Double_t phi)
         :TGeoArb8("", 0, 0)
{
// Constructor providing just a range in Z, theta and phi.
   fDz = dz;
   fTheta = theta;
   fPhi = phi;
   fH1 = fH2 = fBl1 = fBl2 = fTl1 = fTl2 = fAlpha1 = fAlpha2 = 0;
}

//_____________________________________________________________________________
TGeoTrap::TGeoTrap(Double_t dz, Double_t theta, Double_t phi, Double_t h1,
              Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2, 
              Double_t tl2, Double_t alpha2)
         :TGeoArb8("", 0, 0)
{
// Normal constructor.
   fDz = dz;
   fTheta = theta;
   fPhi = phi;
   fH1 = h1;
   fH2 = h2;
   fBl1 = bl1;
   fBl2 = bl2;
   fTl1 = tl1;
   fTl2 = tl2;
   fAlpha1 = alpha1;
   fAlpha2 = alpha2;
   Double_t tx = TMath::Tan(theta*TMath::DegToRad())*TMath::Cos(phi*TMath::DegToRad());
   Double_t ty = TMath::Tan(theta*TMath::DegToRad())*TMath::Sin(phi*TMath::DegToRad());
   Double_t ta1 = TMath::Tan(alpha1*TMath::DegToRad());
   Double_t ta2 = TMath::Tan(alpha2*TMath::DegToRad());
   fXY[0][0] = -dz*tx-h1*ta1-bl1;    fXY[0][1] = -dz*ty-h1;
   fXY[1][0] = -dz*tx+h1*ta1-tl1;    fXY[1][1] = -dz*ty+h1;
   fXY[2][0] = -dz*tx+h1*ta1+tl1;    fXY[2][1] = -dz*ty+h1;
   fXY[3][0] = -dz*tx-h1*ta1+bl1;    fXY[3][1] = -dz*ty-h1;
   fXY[4][0] = dz*tx-h2*ta2-bl2;    fXY[4][1] = dz*ty-h2;
   fXY[5][0] = dz*tx+h2*ta2-tl2;    fXY[5][1] = dz*ty+h2;
   fXY[6][0] = dz*tx+h2*ta2+tl2;    fXY[6][1] = dz*ty+h2;
   fXY[7][0] = dz*tx-h2*ta2+bl2;    fXY[7][1] = dz*ty-h2;
   ComputeTwist();
   if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) || 
       (h2<0) || (bl2<0) || (tl2<0)) {
      SetShapeBit(kGeoRunTimeShape);
   } 
   else TGeoArb8::ComputeBBox();
}

//_____________________________________________________________________________
TGeoTrap::TGeoTrap(const char *name, Double_t dz, Double_t theta, Double_t phi, Double_t h1,
              Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2, 
              Double_t tl2, Double_t alpha2)
         :TGeoArb8(name, 0, 0)
{
// Constructor with name.
   SetName(name);
   fDz = dz;
   fTheta = theta;
   fPhi = phi;
   fH1 = h1;
   fH2 = h2;
   fBl1 = bl1;
   fBl2 = bl2;
   fTl1 = tl1;
   fTl2 = tl2;
   fAlpha1 = alpha1;
   fAlpha2 = alpha2;
   for (Int_t i=0; i<8; i++) {
      fXY[i][0] = 0.0;
      fXY[i][1] = 0.0;
   }   
   Double_t tx = TMath::Tan(theta*TMath::DegToRad())*TMath::Cos(phi*TMath::DegToRad());
   Double_t ty = TMath::Tan(theta*TMath::DegToRad())*TMath::Sin(phi*TMath::DegToRad());
   Double_t ta1 = TMath::Tan(alpha1*TMath::DegToRad());
   Double_t ta2 = TMath::Tan(alpha2*TMath::DegToRad());
   fXY[0][0] = -dz*tx-h1*ta1-bl1;    fXY[0][1] = -dz*ty-h1;
   fXY[1][0] = -dz*tx+h1*ta1-tl1;    fXY[1][1] = -dz*ty+h1;
   fXY[2][0] = -dz*tx+h1*ta1+tl1;    fXY[2][1] = -dz*ty+h1;
   fXY[3][0] = -dz*tx-h1*ta1+bl1;    fXY[3][1] = -dz*ty-h1;
   fXY[4][0] = dz*tx-h2*ta2-bl2;    fXY[4][1] = dz*ty-h2;
   fXY[5][0] = dz*tx+h2*ta2-tl2;    fXY[5][1] = dz*ty+h2;
   fXY[6][0] = dz*tx+h2*ta2+tl2;    fXY[6][1] = dz*ty+h2;
   fXY[7][0] = dz*tx-h2*ta2+bl2;    fXY[7][1] = dz*ty-h2;
   ComputeTwist();
   if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) || 
       (h2<0) || (bl2<0) || (tl2<0)) {
      SetShapeBit(kGeoRunTimeShape);
   } 
   else TGeoArb8::ComputeBBox();
}

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

//_____________________________________________________________________________
Double_t TGeoTrap::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 trapezoid
   if (iact<3 && safe) {
   // compute safe distance
      *safe = Safety(point, kTRUE);
      if (iact==0) return TGeoShape::Big();
      if (iact==1 && step<*safe) return TGeoShape::Big();
   }
   // compute distance to get ouside this shape
//   return TGeoArb8::DistFromInside(point, dir, iact, step, safe);

// compute distance to plane ipl :
// ipl=0 : points 0,4,1,5
// ipl=1 : points 1,5,2,6
// ipl=2 : points 2,6,3,7
// ipl=3 : points 3,7,0,4
   Double_t distmin;
   if (dir[2]<0) {
      distmin=(-fDz-point[2])/dir[2];
   } else {
      if (dir[2]>0) distmin =(fDz-point[2])/dir[2];
      else          distmin = TGeoShape::Big();
   }      
   Double_t xa,xb,xc;
   Double_t ya,yb,yc;
   for (Int_t ipl=0;ipl<4;ipl++) {
      Int_t j = (ipl+1)%4;
      xa=fXY[ipl][0];
      ya=fXY[ipl][1];
      xb=fXY[ipl+4][0];
      yb=fXY[ipl+4][1];
      xc=fXY[j][0];
      yc=fXY[j][1];
      Double_t ax,ay,az;
      ax = xb-xa;
      ay = yb-ya;
      az = 2.*fDz;
      Double_t bx,by;
      bx = xc-xa;
      by = yc-ya;
      Double_t ddotn = -dir[0]*az*by + dir[1]*az*bx+dir[2]*(ax*by-ay*bx);
      if (ddotn<=0) continue; // entering
      Double_t saf = -(point[0]-xa)*az*by + (point[1]-ya)*az*bx + (point[2]+fDz)*(ax*by-ay*bx);
      if (saf>=0.0) return 0.0;
      Double_t s = -saf/ddotn;
      if (s<distmin) distmin=s;
   }
   return distmin;   
}   

//_____________________________________________________________________________
Double_t TGeoTrap::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 trapezoid
   if (iact<3 && safe) {
   // compute safe distance
      *safe = Safety(point, kFALSE);
      if (iact==0) return TGeoShape::Big();
      if (iact==1 && step<*safe) 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();
   // compute distance to get ouside this shape
   Bool_t in = kTRUE;
   Double_t pts[8];
   Double_t snxt;
   Double_t xnew, ynew, znew;
   Int_t i,j;
   if (point[2]<-fDz+TGeoShape::Tolerance()) {
      if (dir[2]<=0) return TGeoShape::Big();
      in = kFALSE;
      snxt = -(fDz+point[2])/dir[2];
      xnew = point[0] + snxt*dir[0];
      ynew = point[1] + snxt*dir[1];
      for (i=0;i<4;i++) {
         j = i<<1;
         pts[j] = fXY[i][0];
         pts[j+1] = fXY[i][1];
      }
      if (InsidePolygon(xnew,ynew,pts)) return snxt;   
   } else if (point[2]>fDz-TGeoShape::Tolerance()) {
      if (dir[2]>=0) return TGeoShape::Big();
      in = kFALSE;
      snxt = (fDz-point[2])/dir[2];
      xnew = point[0] + snxt*dir[0];
      ynew = point[1] + snxt*dir[1];
      for (i=0;i<4;i++) {
         j = i<<1;
         pts[j] = fXY[i+4][0];
         pts[j+1] = fXY[i+4][1];
      }
      if (InsidePolygon(xnew,ynew,pts)) return snxt;    
   }   
   snxt = TGeoShape::Big();   
    

   // check lateral faces
   Double_t dz2 =0.5/fDz;
   Double_t xa,xb,xc,xd;
   Double_t ya,yb,yc,yd;
   Double_t ax,ay,az;
   Double_t bx,by;
   Double_t ddotn, saf;
   Double_t safmin = TGeoShape::Big();
   Bool_t exiting = kFALSE;
   for (i=0; i<4; i++) {
      j = (i+1)%4;
      xa=fXY[i][0];
      ya=fXY[i][1];
      xb=fXY[i+4][0];
      yb=fXY[i+4][1];
      xc=fXY[j][0];
      yc=fXY[j][1];
      xd=fXY[4+j][0];
      yd=fXY[4+j][1];
      ax = xb-xa;
      ay = yb-ya;
      az = 2.*fDz;
      bx = xc-xa;
      by = yc-ya;
      ddotn = -dir[0]*az*by + dir[1]*az*bx+dir[2]*(ax*by-ay*bx);
      saf = (point[0]-xa)*az*by - (point[1]-ya)*az*bx - (point[2]+fDz)*(ax*by-ay*bx);
      
      if (saf<=0) {
         // face visible from point outside
         in = kFALSE;
         if (ddotn>=0) return TGeoShape::Big();      
         snxt = saf/ddotn;
         znew = point[2]+snxt*dir[2];
         if (TMath::Abs(znew)<=fDz) {
            xnew = point[0]+snxt*dir[0];
            ynew = point[1]+snxt*dir[1];
            Double_t tx1 =dz2*(xb-xa);
            Double_t ty1 =dz2*(yb-ya);
            Double_t tx2 =dz2*(xd-xc);
            Double_t ty2 =dz2*(yd-yc);
            Double_t dzp =fDz+znew;
            Double_t xs1 =xa+tx1*dzp;
            Double_t ys1 =ya+ty1*dzp;
            Double_t xs2 =xc+tx2*dzp;
            Double_t ys2 =yc+ty2*dzp;
            if (TMath::Abs(xs1-xs2)>TMath::Abs(ys1-ys2)) {
               if ((xnew-xs1)*(xs2-xnew)>=0) return snxt;
            } else {
               if ((ynew-ys1)*(ys2-ynew)>=0) return snxt;
            }      
         }
      } else {
         if (saf<safmin) {
            safmin = saf;
            if (ddotn>=0) exiting = kTRUE;
            else exiting = kFALSE;
         }   
      }   
   }   
   if (!in) return TGeoShape::Big();
   if (exiting) return TGeoShape::Big();
   return 0.0;
}   

//_____________________________________________________________________________
TGeoVolume *TGeoTrap::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, 
                             Double_t start, Double_t step) 
{
//--- Divide this trapezoid shape belonging to volume "voldiv" into ndiv volumes
// called divname, from start position with the given step. Only Z divisions
// are supported. For Z divisions just return the pointer to the volume to be 
// divided. In case a wrong division axis is supplied, returns pointer to 
// volume that was divided.
   TGeoShape *shape;           //--- shape to be created
   TGeoVolume *vol;            //--- division volume to be created
   TGeoVolumeMulti *vmulti;    //--- generic divided volume
   TGeoPatternFinder *finder;  //--- finder to be attached 
   TString opt = "";           //--- option to be attached
   if (iaxis!=3) {
      Error("Divide", "cannot divide trapezoids on other axis than Z");
      return 0;
   }
   Double_t end = start+ndiv*step;
   Double_t points_lo[8];
   Double_t points_hi[8];
   finder = new TGeoPatternTrapZ(voldiv, ndiv, start, end);
   voldiv->SetFinder(finder);
   finder->SetDivIndex(voldiv->GetNdaughters());
   opt = "Z";
   vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
   Double_t txz = ((TGeoPatternTrapZ*)finder)->GetTxz();
   Double_t tyz = ((TGeoPatternTrapZ*)finder)->GetTyz();
   Double_t zmin, zmax, ox,oy,oz;
   for (Int_t idiv=0; idiv<ndiv; idiv++) {
      zmin = start+idiv*step;
      zmax = start+(idiv+1)*step;
      oz = start+idiv*step+step/2;
      ox = oz*txz;
      oy = oz*tyz;
      SetPlaneVertices(zmin, &points_lo[0]);
      SetPlaneVertices(zmax, &points_hi[0]);
      shape = new TGeoTrap(step/2, fTheta, fPhi);
      for (Int_t vert1=0; vert1<4; vert1++)
         ((TGeoArb8*)shape)->SetVertex(vert1, points_lo[2*vert1]-ox, points_lo[2*vert1+1]-oy);
      for (Int_t vert2=0; vert2<4; vert2++)
         ((TGeoArb8*)shape)->SetVertex(vert2+4, points_hi[2*vert2]-ox, points_hi[2*vert2+1]-oy);
      vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
      vmulti->AddVolume(vol);
      voldiv->AddNodeOffset(vol, idiv, oz, opt.Data());
      ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
   }
   return vmulti;
}   

//_____________________________________________________________________________
TGeoShape *TGeoTrap::GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix * /*mat*/) const
{
// In case shape has some negative parameters, these have to be computed
// in order to fit the mother.
   if (!TestShapeBit(kGeoRunTimeShape)) return 0;
   if (mother->IsRunTimeShape()) {
      Error("GetMakeRuntimeShape", "invalid mother");
      return 0;
   }
   Double_t dz, h1, bl1, tl1, h2, bl2, tl2;
   if (fDz<0) dz=((TGeoTrap*)mother)->GetDz();
   else dz=fDz;

   if (fH1<0)  h1 = ((TGeoTrap*)mother)->GetH1();
   else        h1 = fH1;

   if (fH2<0)  h2 = ((TGeoTrap*)mother)->GetH2();
   else        h2 = fH2;

   if (fBl1<0) bl1 = ((TGeoTrap*)mother)->GetBl1();
   else        bl1 = fBl1;

   if (fBl2<0) bl2 = ((TGeoTrap*)mother)->GetBl2();
   else        bl2 = fBl2;

   if (fTl1<0) tl1 = ((TGeoTrap*)mother)->GetTl1();
   else        tl1 = fTl1;

   if (fTl2<0) tl2 = ((TGeoTrap*)mother)->GetTl2();
   else        tl2 = fTl2;

   return (new TGeoTrap(dz, fTheta, fPhi, h1, bl1, tl1, fAlpha1, h2, bl2, tl2, fAlpha2));
}

//_____________________________________________________________________________
Double_t TGeoTrap::Safety(Double_t *point, Bool_t in) const
{
// Computes the closest distance from given point to this shape.
   Double_t safe = TGeoShape::Big();
   Double_t saf[5];
   Double_t norm[3]; // normal to current facette
   Int_t i,j;       // current facette index
   Double_t x0, y0, z0=-fDz, x1, y1, z1=fDz, x2, y2;
   Double_t ax, ay, az=z1-z0, bx, by;
   Double_t fn;
   //---> compute safety for lateral planes
   for (i=0; i<4; i++) {
      if (in) saf[i] = TGeoShape::Big();
      else    saf[i] = 0.;
      x0 = fXY[i][0];
      y0 = fXY[i][1];
      x1 = fXY[i+4][0];
      y1 = fXY[i+4][1];
      ax = x1-x0;
      ay = y1-y0;
      az = z1-z0;
      j  = (i+1)%4;
      x2 = fXY[j][0];
      y2 = fXY[j][1];
      bx = x2-x0;
      by = y2-y0;
      if (TMath::Abs(bx)<TGeoShape::Tolerance() && TMath::Abs(by)<TGeoShape::Tolerance()) {
         x2 = fXY[4+j][0];
         y2 = fXY[4+j][1];
         bx = x2-x1;
         by = y2-y1;
         if (TMath::Abs(bx)<TGeoShape::Tolerance() && TMath::Abs(by)<TGeoShape::Tolerance()) 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;
      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])/fn; // they should be all positive anyway
      } else {
         saf[i] = -saf[i]/fn;   // only negative values are interesting
      }   
   }
   saf[4] = fDz-TMath::Abs(point[2]);
   if (in) {
      safe = saf[0];
      for (j=1;j<5;j++) if (saf[j] <safe) safe = saf[j];
   } else {
      saf[4]=-saf[4];
      safe = saf[0];
      for (j=1;j<5;j++) if (saf[j] >safe) safe = saf[j];
   }
   return safe;
}

//_____________________________________________________________________________
void TGeoTrap::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 << "   dz     = " << fDz << ";" << endl;
   out << "   theta  = " << fTheta << ";" << endl;
   out << "   phi    = " << fPhi << ";" << endl;
   out << "   h1     = " << fH1<< ";" << endl;
   out << "   bl1    = " << fBl1<< ";" << endl;
   out << "   tl1    = " << fTl1<< ";" << endl;
   out << "   alpha1 = " << fAlpha1 << ";" << endl;
   out << "   h2     = " << fH2 << ";" << endl;
   out << "   bl2    = " << fBl2<< ";" << endl;
   out << "   tl2    = " << fTl2<< ";" << endl;
   out << "   alpha2 = " << fAlpha2 << ";" << endl;
   out << "   TGeoShape *" << GetPointerName() << " = new TGeoTrap(\"" << GetName() << "\", dz,theta,phi,h1,bl1,tl1,alpha1,h2,bl2,tl2,alpha2);" << endl;
   TObject::SetBit(TGeoShape::kGeoSavePrimitive);
}

//_____________________________________________________________________________
void TGeoTrap::SetDimensions(Double_t *param)
{
// Set all arb8 params in one step.
// param[0] = dz
// param[1] = theta
// param[2] = phi
// param[3] = h1
// param[4] = bl1
// param[5] = tl1
// param[6] = alpha1
// param[7] = h2
// param[8] = bl2
// param[9] = tl2
// param[10] = alpha2
   fDz      = param[0];
   fTheta = param[1];
   fPhi   = param[2];
   fH1 = param[3];
   fH2 = param[7];
   fBl1 = param[4];
   fBl2 = param[8];
   fTl1 = param[5];
   fTl2 = param[9];
   fAlpha1 = param[6];
   fAlpha2 = param[10];
   Double_t tx = TMath::Tan(fTheta*TMath::DegToRad())*TMath::Cos(fPhi*TMath::DegToRad());
   Double_t ty = TMath::Tan(fTheta*TMath::DegToRad())*TMath::Sin(fPhi*TMath::DegToRad());
   Double_t ta1 = TMath::Tan(fAlpha1*TMath::DegToRad());
   Double_t ta2 = TMath::Tan(fAlpha2*TMath::DegToRad());
   fXY[0][0] = -fDz*tx-fH1*ta1-fBl1;    fXY[0][1] = -fDz*ty-fH1;
   fXY[1][0] = -fDz*tx+fH1*ta1-fTl1;    fXY[1][1] = -fDz*ty+fH1;
   fXY[2][0] = -fDz*tx+fH1*ta1+fTl1;    fXY[2][1] = -fDz*ty+fH1;
   fXY[3][0] = -fDz*tx-fH1*ta1+fBl1;    fXY[3][1] = -fDz*ty-fH1;
   fXY[4][0] = fDz*tx-fH2*ta2-fBl2;    fXY[4][1] = fDz*ty-fH2;
   fXY[5][0] = fDz*tx+fH2*ta2-fTl2;    fXY[5][1] = fDz*ty+fH2;
   fXY[6][0] = fDz*tx+fH2*ta2+fTl2;    fXY[6][1] = fDz*ty+fH2;
   fXY[7][0] = fDz*tx-fH2*ta2+fBl2;    fXY[7][1] = fDz*ty-fH2;
   ComputeTwist();
   if ((fDz<0) || (fH1<0) || (fBl1<0) || (fTl1<0) || 
       (fH2<0) || (fBl2<0) || (fTl2<0)) {
      SetShapeBit(kGeoRunTimeShape);
   } 
   else TGeoArb8::ComputeBBox();
}   

ClassImp(TGeoGtra)

//_____________________________________________________________________________
TGeoGtra::TGeoGtra()
{
// Default ctor
   fTwistAngle = 0;
}

//_____________________________________________________________________________
TGeoGtra::TGeoGtra(Double_t dz, Double_t theta, Double_t phi, Double_t twist, Double_t h1,
              Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2, 
              Double_t tl2, Double_t alpha2)
         :TGeoTrap(dz, theta, phi, h1, bl1, tl1, alpha1, h2, bl2, tl2, alpha2)     
{
// Constructor. 
   fTheta = theta;
   fPhi = phi;
   fH1 = h1;
   fH2 = h2;
   fBl1 = bl1;
   fBl2 = bl2;
   fTl1 = tl1;
   fTl2 = tl2;
   fAlpha1 = alpha1;
   fAlpha2 = alpha2;
   Double_t x, y, dx, dy, dx1, dx2, th, ph, al1, al2;
   al1 = alpha1*TMath::DegToRad();
   al2 = alpha2*TMath::DegToRad();
   th = theta*TMath::DegToRad();
   ph = phi*TMath::DegToRad();
   dx = 2*dz*TMath::Sin(th)*TMath::Cos(ph);
   dy = 2*dz*TMath::Sin(th)*TMath::Sin(ph);
   fDz = dz;
   dx1 = 2*h1*TMath::Tan(al1);
   dx2 = 2*h2*TMath::Tan(al2);

   fTwistAngle = twist;

   Int_t i;
   for (i=0; i<8; i++) {
      fXY[i][0] = 0.0;
      fXY[i][1] = 0.0;
   }   

   fXY[0][0] = -bl1;                fXY[0][1] = -h1;
   fXY[1][0] = -tl1+dx1;            fXY[1][1] = h1;
   fXY[2][0] = tl1+dx1;             fXY[2][1] = h1;
   fXY[3][0] = bl1;                 fXY[3][1] = -h1;
   fXY[4][0] = -bl2+dx;             fXY[4][1] = -h2+dy;
   fXY[5][0] = -tl2+dx+dx2;         fXY[5][1] = h2+dy;
   fXY[6][0] = tl2+dx+dx2;          fXY[6][1] = h2+dy;
   fXY[7][0] = bl2+dx;              fXY[7][1] = -h2+dy;
   for (i=4; i<8; i++) {
      x = fXY[i][0];
      y = fXY[i][1];
      fXY[i][0] = x*TMath::Cos(twist*TMath::DegToRad()) + y*TMath::Sin(twist*TMath::DegToRad());
      fXY[i][1] = -x*TMath::Sin(twist*TMath::DegToRad()) + y*TMath::Cos(twist*TMath::DegToRad());      
   }
   ComputeTwist();
   if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) || 
       (h2<0) || (bl2<0) || (tl2<0)) SetShapeBit(kGeoRunTimeShape);
   else TGeoArb8::ComputeBBox();
}

//_____________________________________________________________________________
TGeoGtra::TGeoGtra(const char *name, Double_t dz, Double_t theta, Double_t phi, Double_t twist, Double_t h1,
              Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2, 
              Double_t tl2, Double_t alpha2)
         :TGeoTrap(name, dz, theta, phi, h1, bl1, tl1, alpha1, h2, bl2, tl2, alpha2)     
{
// Constructor providing the name of the shape. 
   SetName(name);
   fTheta = theta;
   fPhi = phi;
   fH1 = h1;
   fH2 = h2;
   fBl1 = bl1;
   fBl2 = bl2;
   fTl1 = tl1;
   fTl2 = tl2;
   fAlpha1 = alpha1;
   fAlpha2 = alpha2;
   Double_t x, y, dx, dy, dx1, dx2, th, ph, al1, al2;
   al1 = alpha1*TMath::DegToRad();
   al2 = alpha2*TMath::DegToRad();
   th = theta*TMath::DegToRad();
   ph = phi*TMath::DegToRad();
   dx = 2*dz*TMath::Sin(th)*TMath::Cos(ph);
   dy = 2*dz*TMath::Sin(th)*TMath::Sin(ph);
   fDz = dz;
   dx1 = 2*h1*TMath::Tan(al1);
   dx2 = 2*h2*TMath::Tan(al2);

   fTwistAngle = twist;

   Int_t i;
   for (i=0; i<8; i++) {
      fXY[i][0] = 0.0;
      fXY[i][1] = 0.0;
   }   

   fXY[0][0] = -bl1;                fXY[0][1] = -h1;
   fXY[1][0] = -tl1+dx1;            fXY[1][1] = h1;
   fXY[2][0] = tl1+dx1;             fXY[2][1] = h1;
   fXY[3][0] = bl1;                 fXY[3][1] = -h1;
   fXY[4][0] = -bl2+dx;             fXY[4][1] = -h2+dy;
   fXY[5][0] = -tl2+dx+dx2;         fXY[5][1] = h2+dy;
   fXY[6][0] = tl2+dx+dx2;          fXY[6][1] = h2+dy;
   fXY[7][0] = bl2+dx;              fXY[7][1] = -h2+dy;
   for (i=4; i<8; i++) {
      x = fXY[i][0];
      y = fXY[i][1];
      fXY[i][0] = x*TMath::Cos(twist*TMath::DegToRad()) + y*TMath::Sin(twist*TMath::DegToRad());
      fXY[i][1] = -x*TMath::Sin(twist*TMath::DegToRad()) + y*TMath::Cos(twist*TMath::DegToRad());      
   }
   ComputeTwist();
   if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) || 
       (h2<0) || (bl2<0) || (tl2<0)) SetShapeBit(kGeoRunTimeShape);
   else TGeoArb8::ComputeBBox();
}

//_____________________________________________________________________________
TGeoGtra::~TGeoGtra()
{
// Destructor.
}

//_____________________________________________________________________________
Double_t TGeoGtra::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 shape.
   if (iact<3 && safe) {
   // compute safe distance
      *safe = Safety(point, kTRUE);
      if (iact==0) return TGeoShape::Big();
      if (iact==1 && step<*safe) return TGeoShape::Big();
   }
   // compute distance to get ouside this shape
   return TGeoArb8::DistFromInside(point, dir, iact, step, safe);
}   

//_____________________________________________________________________________
Double_t TGeoGtra::DistFromOutside(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 shape.
   if (iact<3 && safe) {
   // compute safe distance
      *safe = Safety(point, kTRUE);
      if (iact==0) return TGeoShape::Big();
      if (iact==1 && step<*safe) return TGeoShape::Big();
   }
   // compute distance to get ouside this shape
   return TGeoArb8::DistFromOutside(point, dir, iact, step, safe);
}   

//_____________________________________________________________________________
TGeoShape *TGeoGtra::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->IsRunTimeShape()) {
      Error("GetMakeRuntimeShape", "invalid mother");
      return 0;
   }
   Double_t dz, h1, bl1, tl1, h2, bl2, tl2;
   if (fDz<0) dz=((TGeoTrap*)mother)->GetDz();
   else dz=fDz;
   if (fH1<0) 
      h1 = ((TGeoTrap*)mother)->GetH1();
   else 
      h1 = fH1;
   if (fH2<0) 
      h2 = ((TGeoTrap*)mother)->GetH2();
   else 
      h2 = fH2;
   if (fBl1<0) 
      bl1 = ((TGeoTrap*)mother)->GetBl1();
   else 
      bl1 = fBl1;
   if (fBl2<0) 
      bl2 = ((TGeoTrap*)mother)->GetBl2();
   else 
      bl2 = fBl2;
   if (fTl1<0) 
      tl1 = ((TGeoTrap*)mother)->GetTl1();
   else 
      tl1 = fTl1;
   if (fTl2<0) 
      tl2 = ((TGeoTrap*)mother)->GetTl2();
   else 
      tl2 = fTl2;
   return (new TGeoGtra(dz, fTheta, fPhi, fTwistAngle ,h1, bl1, tl1, fAlpha1, h2, bl2, tl2, fAlpha2));
}

//_____________________________________________________________________________
void TGeoGtra::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 << "   dz     = " << fDz << ";" << endl;
   out << "   theta  = " << fTheta << ";" << endl;
   out << "   phi    = " << fPhi << ";" << endl;
   out << "   twist  = " << fTwistAngle << ";" << endl;
   out << "   h1     = " << fH1<< ";" << endl;
   out << "   bl1    = " << fBl1<< ";" << endl;
   out << "   tl1    = " << fTl1<< ";" << endl;
   out << "   alpha1 = " << fAlpha1 << ";" << endl;
   out << "   h2     = " << fH2 << ";" << endl;
   out << "   bl2    = " << fBl2<< ";" << endl;
   out << "   tl2    = " << fTl2<< ";" << endl;
   out << "   alpha2 = " << fAlpha2 << ";" << endl;
   out << "   TGeoShape *" << GetPointerName() << " = new TGeoGtra(\"" << GetName() << "\", dz,theta,phi,twist,h1,bl1,tl1,alpha1,h2,bl2,tl2,alpha2);" << endl;
   TObject::SetBit(TGeoShape::kGeoSavePrimitive);
}

//_____________________________________________________________________________
void TGeoGtra::SetDimensions(Double_t *param)
{
// Set all arb8 params in one step.
// param[0] = dz
// param[1] = theta
// param[2] = phi
// param[3] = h1
// param[4] = bl1
// param[5] = tl1
// param[6] = alpha1
// param[7] = h2
// param[8] = bl2
// param[9] = tl2
// param[10] = alpha2
// param[11] = twist
   fDz      = param[0];
   fTheta = param[1];
   fPhi   = param[2];
   fH1 = param[3];
   fH2 = param[7];
   fBl1 = param[4];
   fBl2 = param[8];
   fTl1 = param[5];
   fTl2 = param[9];
   fAlpha1 = param[6];
   fAlpha2 = param[10];
   fTwistAngle = param[11];
   Double_t x, y, dx, dy, dx1, dx2, th, ph, al1, al2;
   al1 = fAlpha1*TMath::DegToRad();
   al2 = fAlpha2*TMath::DegToRad();
   th = fTheta*TMath::DegToRad();
   ph = fPhi*TMath::DegToRad();
   dx = 2*fDz*TMath::Sin(th)*TMath::Cos(ph);
   dy = 2*fDz*TMath::Sin(th)*TMath::Sin(ph);
   dx1 = 2*fH1*TMath::Tan(al1);
   dx2 = 2*fH2*TMath::Tan(al2);
   Int_t i;
   for (i=0; i<8; i++) {
      fXY[i][0] = 0.0;
      fXY[i][1] = 0.0;
   }   

   fXY[0][0] = -fBl1;                fXY[0][1] = -fH1;
   fXY[1][0] = -fTl1+dx1;            fXY[1][1] = fH1;
   fXY[2][0] = fTl1+dx1;             fXY[2][1] = fH1;
   fXY[3][0] = fBl1;                 fXY[3][1] = -fH1;
   fXY[4][0] = -fBl2+dx;             fXY[4][1] = -fH2+dy;
   fXY[5][0] = -fTl2+dx+dx2;         fXY[5][1] = fH2+dy;
   fXY[6][0] = fTl2+dx+dx2;          fXY[6][1] = fH2+dy;
   fXY[7][0] = fBl2+dx;              fXY[7][1] = -fH2+dy;
   for (i=4; i<8; i++) {
      x = fXY[i][0];
      y = fXY[i][1];
      fXY[i][0] = x*TMath::Cos(fTwistAngle*TMath::DegToRad()) + y*TMath::Sin(fTwistAngle*TMath::DegToRad());
      fXY[i][1] = -x*TMath::Sin(fTwistAngle*TMath::DegToRad()) + y*TMath::Cos(fTwistAngle*TMath::DegToRad());      
   }
   ComputeTwist();
   if ((fDz<0) || (fH1<0) || (fBl1<0) || (fTl1<0) || 
       (fH2<0) || (fBl2<0) || (fTl2<0)) SetShapeBit(kGeoRunTimeShape);
   else TGeoArb8::ComputeBBox();
}   
 TGeoArb8.cxx:1
 TGeoArb8.cxx:2
 TGeoArb8.cxx:3
 TGeoArb8.cxx:4
 TGeoArb8.cxx:5
 TGeoArb8.cxx:6
 TGeoArb8.cxx:7
 TGeoArb8.cxx:8
 TGeoArb8.cxx:9
 TGeoArb8.cxx:10
 TGeoArb8.cxx:11
 TGeoArb8.cxx:12
 TGeoArb8.cxx:13
 TGeoArb8.cxx:14
 TGeoArb8.cxx:15
 TGeoArb8.cxx:16
 TGeoArb8.cxx:17
 TGeoArb8.cxx:18
 TGeoArb8.cxx:19
 TGeoArb8.cxx:20
 TGeoArb8.cxx:21
 TGeoArb8.cxx:22
 TGeoArb8.cxx:23
 TGeoArb8.cxx:24
 TGeoArb8.cxx:25
 TGeoArb8.cxx:26
 TGeoArb8.cxx:27
 TGeoArb8.cxx:28
 TGeoArb8.cxx:29
 TGeoArb8.cxx:30
 TGeoArb8.cxx:31
 TGeoArb8.cxx:32
 TGeoArb8.cxx:33
 TGeoArb8.cxx:34
 TGeoArb8.cxx:35
 TGeoArb8.cxx:36
 TGeoArb8.cxx:37
 TGeoArb8.cxx:38
 TGeoArb8.cxx:39
 TGeoArb8.cxx:40
 TGeoArb8.cxx:41
 TGeoArb8.cxx:42
 TGeoArb8.cxx:43
 TGeoArb8.cxx:44
 TGeoArb8.cxx:45
 TGeoArb8.cxx:46
 TGeoArb8.cxx:47
 TGeoArb8.cxx:48
 TGeoArb8.cxx:49
 TGeoArb8.cxx:50
 TGeoArb8.cxx:51
 TGeoArb8.cxx:52
 TGeoArb8.cxx:53
 TGeoArb8.cxx:54
 TGeoArb8.cxx:55
 TGeoArb8.cxx:56
 TGeoArb8.cxx:57
 TGeoArb8.cxx:58
 TGeoArb8.cxx:59
 TGeoArb8.cxx:60
 TGeoArb8.cxx:61
 TGeoArb8.cxx:62
 TGeoArb8.cxx:63
 TGeoArb8.cxx:64
 TGeoArb8.cxx:65
 TGeoArb8.cxx:66
 TGeoArb8.cxx:67
 TGeoArb8.cxx:68
 TGeoArb8.cxx:69
 TGeoArb8.cxx:70
 TGeoArb8.cxx:71
 TGeoArb8.cxx:72
 TGeoArb8.cxx:73
 TGeoArb8.cxx:74
 TGeoArb8.cxx:75
 TGeoArb8.cxx:76
 TGeoArb8.cxx:77
 TGeoArb8.cxx:78
 TGeoArb8.cxx:79
 TGeoArb8.cxx:80
 TGeoArb8.cxx:81
 TGeoArb8.cxx:82
 TGeoArb8.cxx:83
 TGeoArb8.cxx:84
 TGeoArb8.cxx:85
 TGeoArb8.cxx:86
 TGeoArb8.cxx:87
 TGeoArb8.cxx:88
 TGeoArb8.cxx:89
 TGeoArb8.cxx:90
 TGeoArb8.cxx:91
 TGeoArb8.cxx:92
 TGeoArb8.cxx:93
 TGeoArb8.cxx:94
 TGeoArb8.cxx:95
 TGeoArb8.cxx:96
 TGeoArb8.cxx:97
 TGeoArb8.cxx:98
 TGeoArb8.cxx:99
 TGeoArb8.cxx:100
 TGeoArb8.cxx:101
 TGeoArb8.cxx:102
 TGeoArb8.cxx:103
 TGeoArb8.cxx:104
 TGeoArb8.cxx:105
 TGeoArb8.cxx:106
 TGeoArb8.cxx:107
 TGeoArb8.cxx:108
 TGeoArb8.cxx:109
 TGeoArb8.cxx:110
 TGeoArb8.cxx:111
 TGeoArb8.cxx:112
 TGeoArb8.cxx:113
 TGeoArb8.cxx:114
 TGeoArb8.cxx:115
 TGeoArb8.cxx:116
 TGeoArb8.cxx:117
 TGeoArb8.cxx:118
 TGeoArb8.cxx:119
 TGeoArb8.cxx:120
 TGeoArb8.cxx:121
 TGeoArb8.cxx:122
 TGeoArb8.cxx:123
 TGeoArb8.cxx:124
 TGeoArb8.cxx:125
 TGeoArb8.cxx:126
 TGeoArb8.cxx:127
 TGeoArb8.cxx:128
 TGeoArb8.cxx:129
 TGeoArb8.cxx:130
 TGeoArb8.cxx:131
 TGeoArb8.cxx:132
 TGeoArb8.cxx:133
 TGeoArb8.cxx:134
 TGeoArb8.cxx:135
 TGeoArb8.cxx:136
 TGeoArb8.cxx:137
 TGeoArb8.cxx:138
 TGeoArb8.cxx:139
 TGeoArb8.cxx:140
 TGeoArb8.cxx:141
 TGeoArb8.cxx:142
 TGeoArb8.cxx:143
 TGeoArb8.cxx:144
 TGeoArb8.cxx:145
 TGeoArb8.cxx:146
 TGeoArb8.cxx:147
 TGeoArb8.cxx:148
 TGeoArb8.cxx:149
 TGeoArb8.cxx:150
 TGeoArb8.cxx:151
 TGeoArb8.cxx:152
 TGeoArb8.cxx:153
 TGeoArb8.cxx:154
 TGeoArb8.cxx:155
 TGeoArb8.cxx:156
 TGeoArb8.cxx:157
 TGeoArb8.cxx:158
 TGeoArb8.cxx:159
 TGeoArb8.cxx:160
 TGeoArb8.cxx:161
 TGeoArb8.cxx:162
 TGeoArb8.cxx:163
 TGeoArb8.cxx:164
 TGeoArb8.cxx:165
 TGeoArb8.cxx:166
 TGeoArb8.cxx:167
 TGeoArb8.cxx:168
 TGeoArb8.cxx:169
 TGeoArb8.cxx:170
 TGeoArb8.cxx:171
 TGeoArb8.cxx:172
 TGeoArb8.cxx:173
 TGeoArb8.cxx:174
 TGeoArb8.cxx:175
 TGeoArb8.cxx:176
 TGeoArb8.cxx:177
 TGeoArb8.cxx:178
 TGeoArb8.cxx:179
 TGeoArb8.cxx:180
 TGeoArb8.cxx:181
 TGeoArb8.cxx:182
 TGeoArb8.cxx:183
 TGeoArb8.cxx:184
 TGeoArb8.cxx:185
 TGeoArb8.cxx:186
 TGeoArb8.cxx:187
 TGeoArb8.cxx:188
 TGeoArb8.cxx:189
 TGeoArb8.cxx:190
 TGeoArb8.cxx:191
 TGeoArb8.cxx:192
 TGeoArb8.cxx:193
 TGeoArb8.cxx:194
 TGeoArb8.cxx:195
 TGeoArb8.cxx:196
 TGeoArb8.cxx:197
 TGeoArb8.cxx:198
 TGeoArb8.cxx:199
 TGeoArb8.cxx:200
 TGeoArb8.cxx:201
 TGeoArb8.cxx:202
 TGeoArb8.cxx:203
 TGeoArb8.cxx:204
 TGeoArb8.cxx:205
 TGeoArb8.cxx:206
 TGeoArb8.cxx:207
 TGeoArb8.cxx:208
 TGeoArb8.cxx:209
 TGeoArb8.cxx:210
 TGeoArb8.cxx:211
 TGeoArb8.cxx:212
 TGeoArb8.cxx:213
 TGeoArb8.cxx:214
 TGeoArb8.cxx:215
 TGeoArb8.cxx:216
 TGeoArb8.cxx:217
 TGeoArb8.cxx:218
 TGeoArb8.cxx:219
 TGeoArb8.cxx:220
 TGeoArb8.cxx:221
 TGeoArb8.cxx:222
 TGeoArb8.cxx:223
 TGeoArb8.cxx:224
 TGeoArb8.cxx:225
 TGeoArb8.cxx:226
 TGeoArb8.cxx:227
 TGeoArb8.cxx:228
 TGeoArb8.cxx:229
 TGeoArb8.cxx:230
 TGeoArb8.cxx:231
 TGeoArb8.cxx:232
 TGeoArb8.cxx:233
 TGeoArb8.cxx:234
 TGeoArb8.cxx:235
 TGeoArb8.cxx:236
 TGeoArb8.cxx:237
 TGeoArb8.cxx:238
 TGeoArb8.cxx:239
 TGeoArb8.cxx:240
 TGeoArb8.cxx:241
 TGeoArb8.cxx:242
 TGeoArb8.cxx:243
 TGeoArb8.cxx:244
 TGeoArb8.cxx:245
 TGeoArb8.cxx:246
 TGeoArb8.cxx:247
 TGeoArb8.cxx:248
 TGeoArb8.cxx:249
 TGeoArb8.cxx:250
 TGeoArb8.cxx:251
 TGeoArb8.cxx:252
 TGeoArb8.cxx:253
 TGeoArb8.cxx:254
 TGeoArb8.cxx:255
 TGeoArb8.cxx:256
 TGeoArb8.cxx:257
 TGeoArb8.cxx:258
 TGeoArb8.cxx:259
 TGeoArb8.cxx:260
 TGeoArb8.cxx:261
 TGeoArb8.cxx:262
 TGeoArb8.cxx:263
 TGeoArb8.cxx:264
 TGeoArb8.cxx:265
 TGeoArb8.cxx:266
 TGeoArb8.cxx:267
 TGeoArb8.cxx:268
 TGeoArb8.cxx:269
 TGeoArb8.cxx:270
 TGeoArb8.cxx:271
 TGeoArb8.cxx:272
 TGeoArb8.cxx:273
 TGeoArb8.cxx:274
 TGeoArb8.cxx:275
 TGeoArb8.cxx:276
 TGeoArb8.cxx:277
 TGeoArb8.cxx:278
 TGeoArb8.cxx:279
 TGeoArb8.cxx:280
 TGeoArb8.cxx:281
 TGeoArb8.cxx:282
 TGeoArb8.cxx:283
 TGeoArb8.cxx:284
 TGeoArb8.cxx:285
 TGeoArb8.cxx:286
 TGeoArb8.cxx:287
 TGeoArb8.cxx:288
 TGeoArb8.cxx:289
 TGeoArb8.cxx:290
 TGeoArb8.cxx:291
 TGeoArb8.cxx:292
 TGeoArb8.cxx:293
 TGeoArb8.cxx:294
 TGeoArb8.cxx:295
 TGeoArb8.cxx:296
 TGeoArb8.cxx:297
 TGeoArb8.cxx:298
 TGeoArb8.cxx:299
 TGeoArb8.cxx:300
 TGeoArb8.cxx:301
 TGeoArb8.cxx:302
 TGeoArb8.cxx:303
 TGeoArb8.cxx:304
 TGeoArb8.cxx:305
 TGeoArb8.cxx:306
 TGeoArb8.cxx:307
 TGeoArb8.cxx:308
 TGeoArb8.cxx:309
 TGeoArb8.cxx:310
 TGeoArb8.cxx:311
 TGeoArb8.cxx:312
 TGeoArb8.cxx:313
 TGeoArb8.cxx:314
 TGeoArb8.cxx:315
 TGeoArb8.cxx:316
 TGeoArb8.cxx:317
 TGeoArb8.cxx:318
 TGeoArb8.cxx:319
 TGeoArb8.cxx:320
 TGeoArb8.cxx:321
 TGeoArb8.cxx:322
 TGeoArb8.cxx:323
 TGeoArb8.cxx:324
 TGeoArb8.cxx:325
 TGeoArb8.cxx:326
 TGeoArb8.cxx:327
 TGeoArb8.cxx:328
 TGeoArb8.cxx:329
 TGeoArb8.cxx:330
 TGeoArb8.cxx:331
 TGeoArb8.cxx:332
 TGeoArb8.cxx:333
 TGeoArb8.cxx:334
 TGeoArb8.cxx:335
 TGeoArb8.cxx:336
 TGeoArb8.cxx:337
 TGeoArb8.cxx:338
 TGeoArb8.cxx:339
 TGeoArb8.cxx:340
 TGeoArb8.cxx:341
 TGeoArb8.cxx:342
 TGeoArb8.cxx:343
 TGeoArb8.cxx:344
 TGeoArb8.cxx:345
 TGeoArb8.cxx:346
 TGeoArb8.cxx:347
 TGeoArb8.cxx:348
 TGeoArb8.cxx:349
 TGeoArb8.cxx:350
 TGeoArb8.cxx:351
 TGeoArb8.cxx:352
 TGeoArb8.cxx:353
 TGeoArb8.cxx:354
 TGeoArb8.cxx:355
 TGeoArb8.cxx:356
 TGeoArb8.cxx:357
 TGeoArb8.cxx:358
 TGeoArb8.cxx:359
 TGeoArb8.cxx:360
 TGeoArb8.cxx:361
 TGeoArb8.cxx:362
 TGeoArb8.cxx:363
 TGeoArb8.cxx:364
 TGeoArb8.cxx:365
 TGeoArb8.cxx:366
 TGeoArb8.cxx:367
 TGeoArb8.cxx:368
 TGeoArb8.cxx:369
 TGeoArb8.cxx:370
 TGeoArb8.cxx:371
 TGeoArb8.cxx:372
 TGeoArb8.cxx:373
 TGeoArb8.cxx:374
 TGeoArb8.cxx:375
 TGeoArb8.cxx:376
 TGeoArb8.cxx:377
 TGeoArb8.cxx:378
 TGeoArb8.cxx:379
 TGeoArb8.cxx:380
 TGeoArb8.cxx:381
 TGeoArb8.cxx:382
 TGeoArb8.cxx:383
 TGeoArb8.cxx:384
 TGeoArb8.cxx:385
 TGeoArb8.cxx:386
 TGeoArb8.cxx:387
 TGeoArb8.cxx:388
 TGeoArb8.cxx:389
 TGeoArb8.cxx:390
 TGeoArb8.cxx:391
 TGeoArb8.cxx:392
 TGeoArb8.cxx:393
 TGeoArb8.cxx:394
 TGeoArb8.cxx:395
 TGeoArb8.cxx:396
 TGeoArb8.cxx:397
 TGeoArb8.cxx:398
 TGeoArb8.cxx:399
 TGeoArb8.cxx:400
 TGeoArb8.cxx:401
 TGeoArb8.cxx:402
 TGeoArb8.cxx:403
 TGeoArb8.cxx:404
 TGeoArb8.cxx:405
 TGeoArb8.cxx:406
 TGeoArb8.cxx:407
 TGeoArb8.cxx:408
 TGeoArb8.cxx:409
 TGeoArb8.cxx:410
 TGeoArb8.cxx:411
 TGeoArb8.cxx:412
 TGeoArb8.cxx:413
 TGeoArb8.cxx:414
 TGeoArb8.cxx:415
 TGeoArb8.cxx:416
 TGeoArb8.cxx:417
 TGeoArb8.cxx:418
 TGeoArb8.cxx:419
 TGeoArb8.cxx:420
 TGeoArb8.cxx:421
 TGeoArb8.cxx:422
 TGeoArb8.cxx:423
 TGeoArb8.cxx:424
 TGeoArb8.cxx:425
 TGeoArb8.cxx:426
 TGeoArb8.cxx:427
 TGeoArb8.cxx:428
 TGeoArb8.cxx:429
 TGeoArb8.cxx:430
 TGeoArb8.cxx:431
 TGeoArb8.cxx:432
 TGeoArb8.cxx:433
 TGeoArb8.cxx:434
 TGeoArb8.cxx:435
 TGeoArb8.cxx:436
 TGeoArb8.cxx:437
 TGeoArb8.cxx:438
 TGeoArb8.cxx:439
 TGeoArb8.cxx:440
 TGeoArb8.cxx:441
 TGeoArb8.cxx:442
 TGeoArb8.cxx:443
 TGeoArb8.cxx:444
 TGeoArb8.cxx:445
 TGeoArb8.cxx:446
 TGeoArb8.cxx:447
 TGeoArb8.cxx:448
 TGeoArb8.cxx:449
 TGeoArb8.cxx:450
 TGeoArb8.cxx:451
 TGeoArb8.cxx:452
 TGeoArb8.cxx:453
 TGeoArb8.cxx:454
 TGeoArb8.cxx:455
 TGeoArb8.cxx:456
 TGeoArb8.cxx:457
 TGeoArb8.cxx:458
 TGeoArb8.cxx:459
 TGeoArb8.cxx:460
 TGeoArb8.cxx:461
 TGeoArb8.cxx:462
 TGeoArb8.cxx:463
 TGeoArb8.cxx:464
 TGeoArb8.cxx:465
 TGeoArb8.cxx:466
 TGeoArb8.cxx:467
 TGeoArb8.cxx:468
 TGeoArb8.cxx:469
 TGeoArb8.cxx:470
 TGeoArb8.cxx:471
 TGeoArb8.cxx:472
 TGeoArb8.cxx:473
 TGeoArb8.cxx:474
 TGeoArb8.cxx:475
 TGeoArb8.cxx:476
 TGeoArb8.cxx:477
 TGeoArb8.cxx:478
 TGeoArb8.cxx:479
 TGeoArb8.cxx:480
 TGeoArb8.cxx:481
 TGeoArb8.cxx:482
 TGeoArb8.cxx:483
 TGeoArb8.cxx:484
 TGeoArb8.cxx:485
 TGeoArb8.cxx:486
 TGeoArb8.cxx:487
 TGeoArb8.cxx:488
 TGeoArb8.cxx:489
 TGeoArb8.cxx:490
 TGeoArb8.cxx:491
 TGeoArb8.cxx:492
 TGeoArb8.cxx:493
 TGeoArb8.cxx:494
 TGeoArb8.cxx:495
 TGeoArb8.cxx:496
 TGeoArb8.cxx:497
 TGeoArb8.cxx:498
 TGeoArb8.cxx:499
 TGeoArb8.cxx:500
 TGeoArb8.cxx:501
 TGeoArb8.cxx:502
 TGeoArb8.cxx:503
 TGeoArb8.cxx:504
 TGeoArb8.cxx:505
 TGeoArb8.cxx:506
 TGeoArb8.cxx:507
 TGeoArb8.cxx:508
 TGeoArb8.cxx:509
 TGeoArb8.cxx:510
 TGeoArb8.cxx:511
 TGeoArb8.cxx:512
 TGeoArb8.cxx:513
 TGeoArb8.cxx:514
 TGeoArb8.cxx:515
 TGeoArb8.cxx:516
 TGeoArb8.cxx:517
 TGeoArb8.cxx:518
 TGeoArb8.cxx:519
 TGeoArb8.cxx:520
 TGeoArb8.cxx:521
 TGeoArb8.cxx:522
 TGeoArb8.cxx:523
 TGeoArb8.cxx:524
 TGeoArb8.cxx:525
 TGeoArb8.cxx:526
 TGeoArb8.cxx:527
 TGeoArb8.cxx:528
 TGeoArb8.cxx:529
 TGeoArb8.cxx:530
 TGeoArb8.cxx:531
 TGeoArb8.cxx:532
 TGeoArb8.cxx:533
 TGeoArb8.cxx:534
 TGeoArb8.cxx:535
 TGeoArb8.cxx:536
 TGeoArb8.cxx:537
 TGeoArb8.cxx:538
 TGeoArb8.cxx:539
 TGeoArb8.cxx:540
 TGeoArb8.cxx:541
 TGeoArb8.cxx:542
 TGeoArb8.cxx:543
 TGeoArb8.cxx:544
 TGeoArb8.cxx:545
 TGeoArb8.cxx:546
 TGeoArb8.cxx:547
 TGeoArb8.cxx:548
 TGeoArb8.cxx:549
 TGeoArb8.cxx:550
 TGeoArb8.cxx:551
 TGeoArb8.cxx:552
 TGeoArb8.cxx:553
 TGeoArb8.cxx:554
 TGeoArb8.cxx:555
 TGeoArb8.cxx:556
 TGeoArb8.cxx:557
 TGeoArb8.cxx:558
 TGeoArb8.cxx:559
 TGeoArb8.cxx:560
 TGeoArb8.cxx:561
 TGeoArb8.cxx:562
 TGeoArb8.cxx:563
 TGeoArb8.cxx:564
 TGeoArb8.cxx:565
 TGeoArb8.cxx:566
 TGeoArb8.cxx:567
 TGeoArb8.cxx:568
 TGeoArb8.cxx:569
 TGeoArb8.cxx:570
 TGeoArb8.cxx:571
 TGeoArb8.cxx:572
 TGeoArb8.cxx:573
 TGeoArb8.cxx:574
 TGeoArb8.cxx:575
 TGeoArb8.cxx:576
 TGeoArb8.cxx:577
 TGeoArb8.cxx:578
 TGeoArb8.cxx:579
 TGeoArb8.cxx:580
 TGeoArb8.cxx:581
 TGeoArb8.cxx:582
 TGeoArb8.cxx:583
 TGeoArb8.cxx:584
 TGeoArb8.cxx:585
 TGeoArb8.cxx:586
 TGeoArb8.cxx:587
 TGeoArb8.cxx:588
 TGeoArb8.cxx:589
 TGeoArb8.cxx:590
 TGeoArb8.cxx:591
 TGeoArb8.cxx:592
 TGeoArb8.cxx:593
 TGeoArb8.cxx:594
 TGeoArb8.cxx:595
 TGeoArb8.cxx:596
 TGeoArb8.cxx:597
 TGeoArb8.cxx:598
 TGeoArb8.cxx:599
 TGeoArb8.cxx:600
 TGeoArb8.cxx:601
 TGeoArb8.cxx:602
 TGeoArb8.cxx:603
 TGeoArb8.cxx:604
 TGeoArb8.cxx:605
 TGeoArb8.cxx:606
 TGeoArb8.cxx:607
 TGeoArb8.cxx:608
 TGeoArb8.cxx:609
 TGeoArb8.cxx:610
 TGeoArb8.cxx:611
 TGeoArb8.cxx:612
 TGeoArb8.cxx:613
 TGeoArb8.cxx:614
 TGeoArb8.cxx:615
 TGeoArb8.cxx:616
 TGeoArb8.cxx:617
 TGeoArb8.cxx:618
 TGeoArb8.cxx:619
 TGeoArb8.cxx:620
 TGeoArb8.cxx:621
 TGeoArb8.cxx:622
 TGeoArb8.cxx:623
 TGeoArb8.cxx:624
 TGeoArb8.cxx:625
 TGeoArb8.cxx:626
 TGeoArb8.cxx:627
 TGeoArb8.cxx:628
 TGeoArb8.cxx:629
 TGeoArb8.cxx:630
 TGeoArb8.cxx:631
 TGeoArb8.cxx:632
 TGeoArb8.cxx:633
 TGeoArb8.cxx:634
 TGeoArb8.cxx:635
 TGeoArb8.cxx:636
 TGeoArb8.cxx:637
 TGeoArb8.cxx:638
 TGeoArb8.cxx:639
 TGeoArb8.cxx:640
 TGeoArb8.cxx:641
 TGeoArb8.cxx:642
 TGeoArb8.cxx:643
 TGeoArb8.cxx:644
 TGeoArb8.cxx:645
 TGeoArb8.cxx:646
 TGeoArb8.cxx:647
 TGeoArb8.cxx:648
 TGeoArb8.cxx:649
 TGeoArb8.cxx:650
 TGeoArb8.cxx:651
 TGeoArb8.cxx:652
 TGeoArb8.cxx:653
 TGeoArb8.cxx:654
 TGeoArb8.cxx:655
 TGeoArb8.cxx:656
 TGeoArb8.cxx:657
 TGeoArb8.cxx:658
 TGeoArb8.cxx:659
 TGeoArb8.cxx:660
 TGeoArb8.cxx:661
 TGeoArb8.cxx:662
 TGeoArb8.cxx:663
 TGeoArb8.cxx:664
 TGeoArb8.cxx:665
 TGeoArb8.cxx:666
 TGeoArb8.cxx:667
 TGeoArb8.cxx:668
 TGeoArb8.cxx:669
 TGeoArb8.cxx:670
 TGeoArb8.cxx:671
 TGeoArb8.cxx:672
 TGeoArb8.cxx:673
 TGeoArb8.cxx:674
 TGeoArb8.cxx:675
 TGeoArb8.cxx:676
 TGeoArb8.cxx:677
 TGeoArb8.cxx:678
 TGeoArb8.cxx:679
 TGeoArb8.cxx:680
 TGeoArb8.cxx:681
 TGeoArb8.cxx:682
 TGeoArb8.cxx:683
 TGeoArb8.cxx:684
 TGeoArb8.cxx:685
 TGeoArb8.cxx:686
 TGeoArb8.cxx:687
 TGeoArb8.cxx:688
 TGeoArb8.cxx:689
 TGeoArb8.cxx:690
 TGeoArb8.cxx:691
 TGeoArb8.cxx:692
 TGeoArb8.cxx:693
 TGeoArb8.cxx:694
 TGeoArb8.cxx:695
 TGeoArb8.cxx:696
 TGeoArb8.cxx:697
 TGeoArb8.cxx:698
 TGeoArb8.cxx:699
 TGeoArb8.cxx:700
 TGeoArb8.cxx:701
 TGeoArb8.cxx:702
 TGeoArb8.cxx:703
 TGeoArb8.cxx:704
 TGeoArb8.cxx:705
 TGeoArb8.cxx:706
 TGeoArb8.cxx:707
 TGeoArb8.cxx:708
 TGeoArb8.cxx:709
 TGeoArb8.cxx:710
 TGeoArb8.cxx:711
 TGeoArb8.cxx:712
 TGeoArb8.cxx:713
 TGeoArb8.cxx:714
 TGeoArb8.cxx:715
 TGeoArb8.cxx:716
 TGeoArb8.cxx:717
 TGeoArb8.cxx:718
 TGeoArb8.cxx:719
 TGeoArb8.cxx:720
 TGeoArb8.cxx:721
 TGeoArb8.cxx:722
 TGeoArb8.cxx:723
 TGeoArb8.cxx:724
 TGeoArb8.cxx:725
 TGeoArb8.cxx:726
 TGeoArb8.cxx:727
 TGeoArb8.cxx:728
 TGeoArb8.cxx:729
 TGeoArb8.cxx:730
 TGeoArb8.cxx:731
 TGeoArb8.cxx:732
 TGeoArb8.cxx:733
 TGeoArb8.cxx:734
 TGeoArb8.cxx:735
 TGeoArb8.cxx:736
 TGeoArb8.cxx:737
 TGeoArb8.cxx:738
 TGeoArb8.cxx:739
 TGeoArb8.cxx:740
 TGeoArb8.cxx:741
 TGeoArb8.cxx:742
 TGeoArb8.cxx:743
 TGeoArb8.cxx:744
 TGeoArb8.cxx:745
 TGeoArb8.cxx:746
 TGeoArb8.cxx:747
 TGeoArb8.cxx:748
 TGeoArb8.cxx:749
 TGeoArb8.cxx:750
 TGeoArb8.cxx:751
 TGeoArb8.cxx:752
 TGeoArb8.cxx:753
 TGeoArb8.cxx:754
 TGeoArb8.cxx:755
 TGeoArb8.cxx:756
 TGeoArb8.cxx:757
 TGeoArb8.cxx:758
 TGeoArb8.cxx:759
 TGeoArb8.cxx:760
 TGeoArb8.cxx:761
 TGeoArb8.cxx:762
 TGeoArb8.cxx:763
 TGeoArb8.cxx:764
 TGeoArb8.cxx:765
 TGeoArb8.cxx:766
 TGeoArb8.cxx:767
 TGeoArb8.cxx:768
 TGeoArb8.cxx:769
 TGeoArb8.cxx:770
 TGeoArb8.cxx:771
 TGeoArb8.cxx:772
 TGeoArb8.cxx:773
 TGeoArb8.cxx:774
 TGeoArb8.cxx:775
 TGeoArb8.cxx:776
 TGeoArb8.cxx:777
 TGeoArb8.cxx:778
 TGeoArb8.cxx:779
 TGeoArb8.cxx:780
 TGeoArb8.cxx:781
 TGeoArb8.cxx:782
 TGeoArb8.cxx:783
 TGeoArb8.cxx:784
 TGeoArb8.cxx:785
 TGeoArb8.cxx:786
 TGeoArb8.cxx:787
 TGeoArb8.cxx:788
 TGeoArb8.cxx:789
 TGeoArb8.cxx:790
 TGeoArb8.cxx:791
 TGeoArb8.cxx:792
 TGeoArb8.cxx:793
 TGeoArb8.cxx:794
 TGeoArb8.cxx:795
 TGeoArb8.cxx:796
 TGeoArb8.cxx:797
 TGeoArb8.cxx:798
 TGeoArb8.cxx:799
 TGeoArb8.cxx:800
 TGeoArb8.cxx:801
 TGeoArb8.cxx:802
 TGeoArb8.cxx:803
 TGeoArb8.cxx:804
 TGeoArb8.cxx:805
 TGeoArb8.cxx:806
 TGeoArb8.cxx:807
 TGeoArb8.cxx:808
 TGeoArb8.cxx:809
 TGeoArb8.cxx:810
 TGeoArb8.cxx:811
 TGeoArb8.cxx:812
 TGeoArb8.cxx:813
 TGeoArb8.cxx:814
 TGeoArb8.cxx:815
 TGeoArb8.cxx:816
 TGeoArb8.cxx:817
 TGeoArb8.cxx:818
 TGeoArb8.cxx:819
 TGeoArb8.cxx:820
 TGeoArb8.cxx:821
 TGeoArb8.cxx:822
 TGeoArb8.cxx:823
 TGeoArb8.cxx:824
 TGeoArb8.cxx:825
 TGeoArb8.cxx:826
 TGeoArb8.cxx:827
 TGeoArb8.cxx:828
 TGeoArb8.cxx:829
 TGeoArb8.cxx:830
 TGeoArb8.cxx:831
 TGeoArb8.cxx:832
 TGeoArb8.cxx:833
 TGeoArb8.cxx:834
 TGeoArb8.cxx:835
 TGeoArb8.cxx:836
 TGeoArb8.cxx:837
 TGeoArb8.cxx:838
 TGeoArb8.cxx:839
 TGeoArb8.cxx:840
 TGeoArb8.cxx:841
 TGeoArb8.cxx:842
 TGeoArb8.cxx:843
 TGeoArb8.cxx:844
 TGeoArb8.cxx:845
 TGeoArb8.cxx:846
 TGeoArb8.cxx:847
 TGeoArb8.cxx:848
 TGeoArb8.cxx:849
 TGeoArb8.cxx:850
 TGeoArb8.cxx:851
 TGeoArb8.cxx:852
 TGeoArb8.cxx:853
 TGeoArb8.cxx:854
 TGeoArb8.cxx:855
 TGeoArb8.cxx:856
 TGeoArb8.cxx:857
 TGeoArb8.cxx:858
 TGeoArb8.cxx:859
 TGeoArb8.cxx:860
 TGeoArb8.cxx:861
 TGeoArb8.cxx:862
 TGeoArb8.cxx:863
 TGeoArb8.cxx:864
 TGeoArb8.cxx:865
 TGeoArb8.cxx:866
 TGeoArb8.cxx:867
 TGeoArb8.cxx:868
 TGeoArb8.cxx:869
 TGeoArb8.cxx:870
 TGeoArb8.cxx:871
 TGeoArb8.cxx:872
 TGeoArb8.cxx:873
 TGeoArb8.cxx:874
 TGeoArb8.cxx:875
 TGeoArb8.cxx:876
 TGeoArb8.cxx:877
 TGeoArb8.cxx:878
 TGeoArb8.cxx:879
 TGeoArb8.cxx:880
 TGeoArb8.cxx:881
 TGeoArb8.cxx:882
 TGeoArb8.cxx:883
 TGeoArb8.cxx:884
 TGeoArb8.cxx:885
 TGeoArb8.cxx:886
 TGeoArb8.cxx:887
 TGeoArb8.cxx:888
 TGeoArb8.cxx:889
 TGeoArb8.cxx:890
 TGeoArb8.cxx:891
 TGeoArb8.cxx:892
 TGeoArb8.cxx:893
 TGeoArb8.cxx:894
 TGeoArb8.cxx:895
 TGeoArb8.cxx:896
 TGeoArb8.cxx:897
 TGeoArb8.cxx:898
 TGeoArb8.cxx:899
 TGeoArb8.cxx:900
 TGeoArb8.cxx:901
 TGeoArb8.cxx:902
 TGeoArb8.cxx:903
 TGeoArb8.cxx:904
 TGeoArb8.cxx:905
 TGeoArb8.cxx:906
 TGeoArb8.cxx:907
 TGeoArb8.cxx:908
 TGeoArb8.cxx:909
 TGeoArb8.cxx:910
 TGeoArb8.cxx:911
 TGeoArb8.cxx:912
 TGeoArb8.cxx:913
 TGeoArb8.cxx:914
 TGeoArb8.cxx:915
 TGeoArb8.cxx:916
 TGeoArb8.cxx:917
 TGeoArb8.cxx:918
 TGeoArb8.cxx:919
 TGeoArb8.cxx:920
 TGeoArb8.cxx:921
 TGeoArb8.cxx:922
 TGeoArb8.cxx:923
 TGeoArb8.cxx:924
 TGeoArb8.cxx:925
 TGeoArb8.cxx:926
 TGeoArb8.cxx:927
 TGeoArb8.cxx:928
 TGeoArb8.cxx:929
 TGeoArb8.cxx:930
 TGeoArb8.cxx:931
 TGeoArb8.cxx:932
 TGeoArb8.cxx:933
 TGeoArb8.cxx:934
 TGeoArb8.cxx:935
 TGeoArb8.cxx:936
 TGeoArb8.cxx:937
 TGeoArb8.cxx:938
 TGeoArb8.cxx:939
 TGeoArb8.cxx:940
 TGeoArb8.cxx:941
 TGeoArb8.cxx:942
 TGeoArb8.cxx:943
 TGeoArb8.cxx:944
 TGeoArb8.cxx:945
 TGeoArb8.cxx:946
 TGeoArb8.cxx:947
 TGeoArb8.cxx:948
 TGeoArb8.cxx:949
 TGeoArb8.cxx:950
 TGeoArb8.cxx:951
 TGeoArb8.cxx:952
 TGeoArb8.cxx:953
 TGeoArb8.cxx:954
 TGeoArb8.cxx:955
 TGeoArb8.cxx:956
 TGeoArb8.cxx:957
 TGeoArb8.cxx:958
 TGeoArb8.cxx:959
 TGeoArb8.cxx:960
 TGeoArb8.cxx:961
 TGeoArb8.cxx:962
 TGeoArb8.cxx:963
 TGeoArb8.cxx:964
 TGeoArb8.cxx:965
 TGeoArb8.cxx:966
 TGeoArb8.cxx:967
 TGeoArb8.cxx:968
 TGeoArb8.cxx:969
 TGeoArb8.cxx:970
 TGeoArb8.cxx:971
 TGeoArb8.cxx:972
 TGeoArb8.cxx:973
 TGeoArb8.cxx:974
 TGeoArb8.cxx:975
 TGeoArb8.cxx:976
 TGeoArb8.cxx:977
 TGeoArb8.cxx:978
 TGeoArb8.cxx:979
 TGeoArb8.cxx:980
 TGeoArb8.cxx:981
 TGeoArb8.cxx:982
 TGeoArb8.cxx:983
 TGeoArb8.cxx:984
 TGeoArb8.cxx:985
 TGeoArb8.cxx:986
 TGeoArb8.cxx:987
 TGeoArb8.cxx:988
 TGeoArb8.cxx:989
 TGeoArb8.cxx:990
 TGeoArb8.cxx:991
 TGeoArb8.cxx:992
 TGeoArb8.cxx:993
 TGeoArb8.cxx:994
 TGeoArb8.cxx:995
 TGeoArb8.cxx:996
 TGeoArb8.cxx:997
 TGeoArb8.cxx:998
 TGeoArb8.cxx:999
 TGeoArb8.cxx:1000
 TGeoArb8.cxx:1001
 TGeoArb8.cxx:1002
 TGeoArb8.cxx:1003
 TGeoArb8.cxx:1004
 TGeoArb8.cxx:1005
 TGeoArb8.cxx:1006
 TGeoArb8.cxx:1007
 TGeoArb8.cxx:1008
 TGeoArb8.cxx:1009
 TGeoArb8.cxx:1010
 TGeoArb8.cxx:1011
 TGeoArb8.cxx:1012
 TGeoArb8.cxx:1013
 TGeoArb8.cxx:1014
 TGeoArb8.cxx:1015
 TGeoArb8.cxx:1016
 TGeoArb8.cxx:1017
 TGeoArb8.cxx:1018
 TGeoArb8.cxx:1019
 TGeoArb8.cxx:1020
 TGeoArb8.cxx:1021
 TGeoArb8.cxx:1022
 TGeoArb8.cxx:1023
 TGeoArb8.cxx:1024
 TGeoArb8.cxx:1025
 TGeoArb8.cxx:1026
 TGeoArb8.cxx:1027
 TGeoArb8.cxx:1028
 TGeoArb8.cxx:1029
 TGeoArb8.cxx:1030
 TGeoArb8.cxx:1031
 TGeoArb8.cxx:1032
 TGeoArb8.cxx:1033
 TGeoArb8.cxx:1034
 TGeoArb8.cxx:1035
 TGeoArb8.cxx:1036
 TGeoArb8.cxx:1037
 TGeoArb8.cxx:1038
 TGeoArb8.cxx:1039
 TGeoArb8.cxx:1040
 TGeoArb8.cxx:1041
 TGeoArb8.cxx:1042
 TGeoArb8.cxx:1043
 TGeoArb8.cxx:1044
 TGeoArb8.cxx:1045
 TGeoArb8.cxx:1046
 TGeoArb8.cxx:1047
 TGeoArb8.cxx:1048
 TGeoArb8.cxx:1049
 TGeoArb8.cxx:1050
 TGeoArb8.cxx:1051
 TGeoArb8.cxx:1052
 TGeoArb8.cxx:1053
 TGeoArb8.cxx:1054
 TGeoArb8.cxx:1055
 TGeoArb8.cxx:1056
 TGeoArb8.cxx:1057
 TGeoArb8.cxx:1058
 TGeoArb8.cxx:1059
 TGeoArb8.cxx:1060
 TGeoArb8.cxx:1061
 TGeoArb8.cxx:1062
 TGeoArb8.cxx:1063
 TGeoArb8.cxx:1064
 TGeoArb8.cxx:1065
 TGeoArb8.cxx:1066
 TGeoArb8.cxx:1067
 TGeoArb8.cxx:1068
 TGeoArb8.cxx:1069
 TGeoArb8.cxx:1070
 TGeoArb8.cxx:1071
 TGeoArb8.cxx:1072
 TGeoArb8.cxx:1073
 TGeoArb8.cxx:1074
 TGeoArb8.cxx:1075
 TGeoArb8.cxx:1076
 TGeoArb8.cxx:1077
 TGeoArb8.cxx:1078
 TGeoArb8.cxx:1079
 TGeoArb8.cxx:1080
 TGeoArb8.cxx:1081
 TGeoArb8.cxx:1082
 TGeoArb8.cxx:1083
 TGeoArb8.cxx:1084
 TGeoArb8.cxx:1085
 TGeoArb8.cxx:1086
 TGeoArb8.cxx:1087
 TGeoArb8.cxx:1088
 TGeoArb8.cxx:1089
 TGeoArb8.cxx:1090
 TGeoArb8.cxx:1091
 TGeoArb8.cxx:1092
 TGeoArb8.cxx:1093
 TGeoArb8.cxx:1094
 TGeoArb8.cxx:1095
 TGeoArb8.cxx:1096
 TGeoArb8.cxx:1097
 TGeoArb8.cxx:1098
 TGeoArb8.cxx:1099
 TGeoArb8.cxx:1100
 TGeoArb8.cxx:1101
 TGeoArb8.cxx:1102
 TGeoArb8.cxx:1103
 TGeoArb8.cxx:1104
 TGeoArb8.cxx:1105
 TGeoArb8.cxx:1106
 TGeoArb8.cxx:1107
 TGeoArb8.cxx:1108
 TGeoArb8.cxx:1109
 TGeoArb8.cxx:1110
 TGeoArb8.cxx:1111
 TGeoArb8.cxx:1112
 TGeoArb8.cxx:1113
 TGeoArb8.cxx:1114
 TGeoArb8.cxx:1115
 TGeoArb8.cxx:1116
 TGeoArb8.cxx:1117
 TGeoArb8.cxx:1118
 TGeoArb8.cxx:1119
 TGeoArb8.cxx:1120
 TGeoArb8.cxx:1121
 TGeoArb8.cxx:1122
 TGeoArb8.cxx:1123
 TGeoArb8.cxx:1124
 TGeoArb8.cxx:1125
 TGeoArb8.cxx:1126
 TGeoArb8.cxx:1127
 TGeoArb8.cxx:1128
 TGeoArb8.cxx:1129
 TGeoArb8.cxx:1130
 TGeoArb8.cxx:1131
 TGeoArb8.cxx:1132
 TGeoArb8.cxx:1133
 TGeoArb8.cxx:1134
 TGeoArb8.cxx:1135
 TGeoArb8.cxx:1136
 TGeoArb8.cxx:1137
 TGeoArb8.cxx:1138
 TGeoArb8.cxx:1139
 TGeoArb8.cxx:1140
 TGeoArb8.cxx:1141
 TGeoArb8.cxx:1142
 TGeoArb8.cxx:1143
 TGeoArb8.cxx:1144
 TGeoArb8.cxx:1145
 TGeoArb8.cxx:1146
 TGeoArb8.cxx:1147
 TGeoArb8.cxx:1148
 TGeoArb8.cxx:1149
 TGeoArb8.cxx:1150
 TGeoArb8.cxx:1151
 TGeoArb8.cxx:1152
 TGeoArb8.cxx:1153
 TGeoArb8.cxx:1154
 TGeoArb8.cxx:1155
 TGeoArb8.cxx:1156
 TGeoArb8.cxx:1157
 TGeoArb8.cxx:1158
 TGeoArb8.cxx:1159
 TGeoArb8.cxx:1160
 TGeoArb8.cxx:1161
 TGeoArb8.cxx:1162
 TGeoArb8.cxx:1163
 TGeoArb8.cxx:1164
 TGeoArb8.cxx:1165
 TGeoArb8.cxx:1166
 TGeoArb8.cxx:1167
 TGeoArb8.cxx:1168
 TGeoArb8.cxx:1169
 TGeoArb8.cxx:1170
 TGeoArb8.cxx:1171
 TGeoArb8.cxx:1172
 TGeoArb8.cxx:1173
 TGeoArb8.cxx:1174
 TGeoArb8.cxx:1175
 TGeoArb8.cxx:1176
 TGeoArb8.cxx:1177
 TGeoArb8.cxx:1178
 TGeoArb8.cxx:1179
 TGeoArb8.cxx:1180
 TGeoArb8.cxx:1181
 TGeoArb8.cxx:1182
 TGeoArb8.cxx:1183
 TGeoArb8.cxx:1184
 TGeoArb8.cxx:1185
 TGeoArb8.cxx:1186
 TGeoArb8.cxx:1187
 TGeoArb8.cxx:1188
 TGeoArb8.cxx:1189
 TGeoArb8.cxx:1190
 TGeoArb8.cxx:1191
 TGeoArb8.cxx:1192
 TGeoArb8.cxx:1193
 TGeoArb8.cxx:1194
 TGeoArb8.cxx:1195
 TGeoArb8.cxx:1196
 TGeoArb8.cxx:1197
 TGeoArb8.cxx:1198
 TGeoArb8.cxx:1199
 TGeoArb8.cxx:1200
 TGeoArb8.cxx:1201
 TGeoArb8.cxx:1202
 TGeoArb8.cxx:1203
 TGeoArb8.cxx:1204
 TGeoArb8.cxx:1205
 TGeoArb8.cxx:1206
 TGeoArb8.cxx:1207
 TGeoArb8.cxx:1208
 TGeoArb8.cxx:1209
 TGeoArb8.cxx:1210
 TGeoArb8.cxx:1211
 TGeoArb8.cxx:1212
 TGeoArb8.cxx:1213
 TGeoArb8.cxx:1214
 TGeoArb8.cxx:1215
 TGeoArb8.cxx:1216
 TGeoArb8.cxx:1217
 TGeoArb8.cxx:1218
 TGeoArb8.cxx:1219
 TGeoArb8.cxx:1220
 TGeoArb8.cxx:1221
 TGeoArb8.cxx:1222
 TGeoArb8.cxx:1223
 TGeoArb8.cxx:1224
 TGeoArb8.cxx:1225
 TGeoArb8.cxx:1226
 TGeoArb8.cxx:1227
 TGeoArb8.cxx:1228
 TGeoArb8.cxx:1229
 TGeoArb8.cxx:1230
 TGeoArb8.cxx:1231
 TGeoArb8.cxx:1232
 TGeoArb8.cxx:1233
 TGeoArb8.cxx:1234
 TGeoArb8.cxx:1235
 TGeoArb8.cxx:1236
 TGeoArb8.cxx:1237
 TGeoArb8.cxx:1238
 TGeoArb8.cxx:1239
 TGeoArb8.cxx:1240
 TGeoArb8.cxx:1241
 TGeoArb8.cxx:1242
 TGeoArb8.cxx:1243
 TGeoArb8.cxx:1244
 TGeoArb8.cxx:1245
 TGeoArb8.cxx:1246
 TGeoArb8.cxx:1247
 TGeoArb8.cxx:1248
 TGeoArb8.cxx:1249
 TGeoArb8.cxx:1250
 TGeoArb8.cxx:1251
 TGeoArb8.cxx:1252
 TGeoArb8.cxx:1253
 TGeoArb8.cxx:1254
 TGeoArb8.cxx:1255
 TGeoArb8.cxx:1256
 TGeoArb8.cxx:1257
 TGeoArb8.cxx:1258
 TGeoArb8.cxx:1259
 TGeoArb8.cxx:1260
 TGeoArb8.cxx:1261
 TGeoArb8.cxx:1262
 TGeoArb8.cxx:1263
 TGeoArb8.cxx:1264
 TGeoArb8.cxx:1265
 TGeoArb8.cxx:1266
 TGeoArb8.cxx:1267
 TGeoArb8.cxx:1268
 TGeoArb8.cxx:1269
 TGeoArb8.cxx:1270
 TGeoArb8.cxx:1271
 TGeoArb8.cxx:1272
 TGeoArb8.cxx:1273
 TGeoArb8.cxx:1274
 TGeoArb8.cxx:1275
 TGeoArb8.cxx:1276
 TGeoArb8.cxx:1277
 TGeoArb8.cxx:1278
 TGeoArb8.cxx:1279
 TGeoArb8.cxx:1280
 TGeoArb8.cxx:1281
 TGeoArb8.cxx:1282
 TGeoArb8.cxx:1283
 TGeoArb8.cxx:1284
 TGeoArb8.cxx:1285
 TGeoArb8.cxx:1286
 TGeoArb8.cxx:1287
 TGeoArb8.cxx:1288
 TGeoArb8.cxx:1289
 TGeoArb8.cxx:1290
 TGeoArb8.cxx:1291
 TGeoArb8.cxx:1292
 TGeoArb8.cxx:1293
 TGeoArb8.cxx:1294
 TGeoArb8.cxx:1295
 TGeoArb8.cxx:1296
 TGeoArb8.cxx:1297
 TGeoArb8.cxx:1298
 TGeoArb8.cxx:1299
 TGeoArb8.cxx:1300
 TGeoArb8.cxx:1301
 TGeoArb8.cxx:1302
 TGeoArb8.cxx:1303
 TGeoArb8.cxx:1304
 TGeoArb8.cxx:1305
 TGeoArb8.cxx:1306
 TGeoArb8.cxx:1307
 TGeoArb8.cxx:1308
 TGeoArb8.cxx:1309
 TGeoArb8.cxx:1310
 TGeoArb8.cxx:1311
 TGeoArb8.cxx:1312
 TGeoArb8.cxx:1313
 TGeoArb8.cxx:1314
 TGeoArb8.cxx:1315
 TGeoArb8.cxx:1316
 TGeoArb8.cxx:1317
 TGeoArb8.cxx:1318
 TGeoArb8.cxx:1319
 TGeoArb8.cxx:1320
 TGeoArb8.cxx:1321
 TGeoArb8.cxx:1322
 TGeoArb8.cxx:1323
 TGeoArb8.cxx:1324
 TGeoArb8.cxx:1325
 TGeoArb8.cxx:1326
 TGeoArb8.cxx:1327
 TGeoArb8.cxx:1328
 TGeoArb8.cxx:1329
 TGeoArb8.cxx:1330
 TGeoArb8.cxx:1331
 TGeoArb8.cxx:1332
 TGeoArb8.cxx:1333
 TGeoArb8.cxx:1334
 TGeoArb8.cxx:1335
 TGeoArb8.cxx:1336
 TGeoArb8.cxx:1337
 TGeoArb8.cxx:1338
 TGeoArb8.cxx:1339
 TGeoArb8.cxx:1340
 TGeoArb8.cxx:1341
 TGeoArb8.cxx:1342
 TGeoArb8.cxx:1343
 TGeoArb8.cxx:1344
 TGeoArb8.cxx:1345
 TGeoArb8.cxx:1346
 TGeoArb8.cxx:1347
 TGeoArb8.cxx:1348
 TGeoArb8.cxx:1349
 TGeoArb8.cxx:1350
 TGeoArb8.cxx:1351
 TGeoArb8.cxx:1352
 TGeoArb8.cxx:1353
 TGeoArb8.cxx:1354
 TGeoArb8.cxx:1355
 TGeoArb8.cxx:1356
 TGeoArb8.cxx:1357
 TGeoArb8.cxx:1358
 TGeoArb8.cxx:1359
 TGeoArb8.cxx:1360
 TGeoArb8.cxx:1361
 TGeoArb8.cxx:1362
 TGeoArb8.cxx:1363
 TGeoArb8.cxx:1364
 TGeoArb8.cxx:1365
 TGeoArb8.cxx:1366
 TGeoArb8.cxx:1367
 TGeoArb8.cxx:1368
 TGeoArb8.cxx:1369
 TGeoArb8.cxx:1370
 TGeoArb8.cxx:1371
 TGeoArb8.cxx:1372
 TGeoArb8.cxx:1373
 TGeoArb8.cxx:1374
 TGeoArb8.cxx:1375
 TGeoArb8.cxx:1376
 TGeoArb8.cxx:1377
 TGeoArb8.cxx:1378
 TGeoArb8.cxx:1379
 TGeoArb8.cxx:1380
 TGeoArb8.cxx:1381
 TGeoArb8.cxx:1382
 TGeoArb8.cxx:1383
 TGeoArb8.cxx:1384
 TGeoArb8.cxx:1385
 TGeoArb8.cxx:1386
 TGeoArb8.cxx:1387
 TGeoArb8.cxx:1388
 TGeoArb8.cxx:1389
 TGeoArb8.cxx:1390
 TGeoArb8.cxx:1391
 TGeoArb8.cxx:1392
 TGeoArb8.cxx:1393
 TGeoArb8.cxx:1394
 TGeoArb8.cxx:1395
 TGeoArb8.cxx:1396
 TGeoArb8.cxx:1397
 TGeoArb8.cxx:1398
 TGeoArb8.cxx:1399
 TGeoArb8.cxx:1400
 TGeoArb8.cxx:1401
 TGeoArb8.cxx:1402
 TGeoArb8.cxx:1403
 TGeoArb8.cxx:1404
 TGeoArb8.cxx:1405
 TGeoArb8.cxx:1406
 TGeoArb8.cxx:1407
 TGeoArb8.cxx:1408
 TGeoArb8.cxx:1409
 TGeoArb8.cxx:1410
 TGeoArb8.cxx:1411
 TGeoArb8.cxx:1412
 TGeoArb8.cxx:1413
 TGeoArb8.cxx:1414
 TGeoArb8.cxx:1415
 TGeoArb8.cxx:1416
 TGeoArb8.cxx:1417
 TGeoArb8.cxx:1418
 TGeoArb8.cxx:1419
 TGeoArb8.cxx:1420
 TGeoArb8.cxx:1421
 TGeoArb8.cxx:1422
 TGeoArb8.cxx:1423
 TGeoArb8.cxx:1424
 TGeoArb8.cxx:1425
 TGeoArb8.cxx:1426
 TGeoArb8.cxx:1427
 TGeoArb8.cxx:1428
 TGeoArb8.cxx:1429
 TGeoArb8.cxx:1430
 TGeoArb8.cxx:1431
 TGeoArb8.cxx:1432
 TGeoArb8.cxx:1433
 TGeoArb8.cxx:1434
 TGeoArb8.cxx:1435
 TGeoArb8.cxx:1436
 TGeoArb8.cxx:1437
 TGeoArb8.cxx:1438
 TGeoArb8.cxx:1439
 TGeoArb8.cxx:1440
 TGeoArb8.cxx:1441
 TGeoArb8.cxx:1442
 TGeoArb8.cxx:1443
 TGeoArb8.cxx:1444
 TGeoArb8.cxx:1445
 TGeoArb8.cxx:1446
 TGeoArb8.cxx:1447
 TGeoArb8.cxx:1448
 TGeoArb8.cxx:1449
 TGeoArb8.cxx:1450
 TGeoArb8.cxx:1451
 TGeoArb8.cxx:1452
 TGeoArb8.cxx:1453
 TGeoArb8.cxx:1454
 TGeoArb8.cxx:1455
 TGeoArb8.cxx:1456
 TGeoArb8.cxx:1457
 TGeoArb8.cxx:1458
 TGeoArb8.cxx:1459
 TGeoArb8.cxx:1460
 TGeoArb8.cxx:1461
 TGeoArb8.cxx:1462
 TGeoArb8.cxx:1463
 TGeoArb8.cxx:1464
 TGeoArb8.cxx:1465
 TGeoArb8.cxx:1466
 TGeoArb8.cxx:1467
 TGeoArb8.cxx:1468
 TGeoArb8.cxx:1469
 TGeoArb8.cxx:1470
 TGeoArb8.cxx:1471
 TGeoArb8.cxx:1472
 TGeoArb8.cxx:1473
 TGeoArb8.cxx:1474
 TGeoArb8.cxx:1475
 TGeoArb8.cxx:1476
 TGeoArb8.cxx:1477
 TGeoArb8.cxx:1478
 TGeoArb8.cxx:1479
 TGeoArb8.cxx:1480
 TGeoArb8.cxx:1481
 TGeoArb8.cxx:1482
 TGeoArb8.cxx:1483
 TGeoArb8.cxx:1484
 TGeoArb8.cxx:1485
 TGeoArb8.cxx:1486
 TGeoArb8.cxx:1487
 TGeoArb8.cxx:1488
 TGeoArb8.cxx:1489
 TGeoArb8.cxx:1490
 TGeoArb8.cxx:1491
 TGeoArb8.cxx:1492
 TGeoArb8.cxx:1493
 TGeoArb8.cxx:1494
 TGeoArb8.cxx:1495
 TGeoArb8.cxx:1496
 TGeoArb8.cxx:1497
 TGeoArb8.cxx:1498
 TGeoArb8.cxx:1499
 TGeoArb8.cxx:1500
 TGeoArb8.cxx:1501
 TGeoArb8.cxx:1502
 TGeoArb8.cxx:1503
 TGeoArb8.cxx:1504
 TGeoArb8.cxx:1505
 TGeoArb8.cxx:1506
 TGeoArb8.cxx:1507
 TGeoArb8.cxx:1508
 TGeoArb8.cxx:1509
 TGeoArb8.cxx:1510
 TGeoArb8.cxx:1511
 TGeoArb8.cxx:1512
 TGeoArb8.cxx:1513
 TGeoArb8.cxx:1514
 TGeoArb8.cxx:1515
 TGeoArb8.cxx:1516
 TGeoArb8.cxx:1517
 TGeoArb8.cxx:1518
 TGeoArb8.cxx:1519
 TGeoArb8.cxx:1520
 TGeoArb8.cxx:1521
 TGeoArb8.cxx:1522
 TGeoArb8.cxx:1523
 TGeoArb8.cxx:1524
 TGeoArb8.cxx:1525
 TGeoArb8.cxx:1526
 TGeoArb8.cxx:1527
 TGeoArb8.cxx:1528
 TGeoArb8.cxx:1529
 TGeoArb8.cxx:1530
 TGeoArb8.cxx:1531
 TGeoArb8.cxx:1532
 TGeoArb8.cxx:1533
 TGeoArb8.cxx:1534
 TGeoArb8.cxx:1535
 TGeoArb8.cxx:1536
 TGeoArb8.cxx:1537
 TGeoArb8.cxx:1538
 TGeoArb8.cxx:1539
 TGeoArb8.cxx:1540
 TGeoArb8.cxx:1541
 TGeoArb8.cxx:1542
 TGeoArb8.cxx:1543
 TGeoArb8.cxx:1544
 TGeoArb8.cxx:1545
 TGeoArb8.cxx:1546
 TGeoArb8.cxx:1547
 TGeoArb8.cxx:1548
 TGeoArb8.cxx:1549
 TGeoArb8.cxx:1550
 TGeoArb8.cxx:1551
 TGeoArb8.cxx:1552
 TGeoArb8.cxx:1553
 TGeoArb8.cxx:1554
 TGeoArb8.cxx:1555
 TGeoArb8.cxx:1556
 TGeoArb8.cxx:1557
 TGeoArb8.cxx:1558
 TGeoArb8.cxx:1559
 TGeoArb8.cxx:1560
 TGeoArb8.cxx:1561
 TGeoArb8.cxx:1562
 TGeoArb8.cxx:1563
 TGeoArb8.cxx:1564
 TGeoArb8.cxx:1565
 TGeoArb8.cxx:1566
 TGeoArb8.cxx:1567
 TGeoArb8.cxx:1568
 TGeoArb8.cxx:1569
 TGeoArb8.cxx:1570
 TGeoArb8.cxx:1571
 TGeoArb8.cxx:1572
 TGeoArb8.cxx:1573
 TGeoArb8.cxx:1574
 TGeoArb8.cxx:1575
 TGeoArb8.cxx:1576
 TGeoArb8.cxx:1577
 TGeoArb8.cxx:1578
 TGeoArb8.cxx:1579
 TGeoArb8.cxx:1580
 TGeoArb8.cxx:1581
 TGeoArb8.cxx:1582
 TGeoArb8.cxx:1583
 TGeoArb8.cxx:1584
 TGeoArb8.cxx:1585
 TGeoArb8.cxx:1586
 TGeoArb8.cxx:1587
 TGeoArb8.cxx:1588
 TGeoArb8.cxx:1589
 TGeoArb8.cxx:1590
 TGeoArb8.cxx:1591
 TGeoArb8.cxx:1592
 TGeoArb8.cxx:1593
 TGeoArb8.cxx:1594
 TGeoArb8.cxx:1595
 TGeoArb8.cxx:1596
 TGeoArb8.cxx:1597
 TGeoArb8.cxx:1598
 TGeoArb8.cxx:1599
 TGeoArb8.cxx:1600
 TGeoArb8.cxx:1601
 TGeoArb8.cxx:1602
 TGeoArb8.cxx:1603
 TGeoArb8.cxx:1604
 TGeoArb8.cxx:1605
 TGeoArb8.cxx:1606
 TGeoArb8.cxx:1607
 TGeoArb8.cxx:1608
 TGeoArb8.cxx:1609
 TGeoArb8.cxx:1610
 TGeoArb8.cxx:1611
 TGeoArb8.cxx:1612
 TGeoArb8.cxx:1613
 TGeoArb8.cxx:1614
 TGeoArb8.cxx:1615
 TGeoArb8.cxx:1616
 TGeoArb8.cxx:1617
 TGeoArb8.cxx:1618
 TGeoArb8.cxx:1619
 TGeoArb8.cxx:1620
 TGeoArb8.cxx:1621
 TGeoArb8.cxx:1622
 TGeoArb8.cxx:1623
 TGeoArb8.cxx:1624
 TGeoArb8.cxx:1625
 TGeoArb8.cxx:1626
 TGeoArb8.cxx:1627
 TGeoArb8.cxx:1628
 TGeoArb8.cxx:1629
 TGeoArb8.cxx:1630
 TGeoArb8.cxx:1631
 TGeoArb8.cxx:1632
 TGeoArb8.cxx:1633
 TGeoArb8.cxx:1634
 TGeoArb8.cxx:1635
 TGeoArb8.cxx:1636
 TGeoArb8.cxx:1637
 TGeoArb8.cxx:1638
 TGeoArb8.cxx:1639
 TGeoArb8.cxx:1640
 TGeoArb8.cxx:1641
 TGeoArb8.cxx:1642
 TGeoArb8.cxx:1643
 TGeoArb8.cxx:1644
 TGeoArb8.cxx:1645
 TGeoArb8.cxx:1646
 TGeoArb8.cxx:1647
 TGeoArb8.cxx:1648
 TGeoArb8.cxx:1649
 TGeoArb8.cxx:1650
 TGeoArb8.cxx:1651
 TGeoArb8.cxx:1652
 TGeoArb8.cxx:1653
 TGeoArb8.cxx:1654
 TGeoArb8.cxx:1655
 TGeoArb8.cxx:1656
 TGeoArb8.cxx:1657
 TGeoArb8.cxx:1658
 TGeoArb8.cxx:1659
 TGeoArb8.cxx:1660
 TGeoArb8.cxx:1661
 TGeoArb8.cxx:1662
 TGeoArb8.cxx:1663
 TGeoArb8.cxx:1664
 TGeoArb8.cxx:1665
 TGeoArb8.cxx:1666
 TGeoArb8.cxx:1667
 TGeoArb8.cxx:1668
 TGeoArb8.cxx:1669
 TGeoArb8.cxx:1670
 TGeoArb8.cxx:1671
 TGeoArb8.cxx:1672
 TGeoArb8.cxx:1673
 TGeoArb8.cxx:1674
 TGeoArb8.cxx:1675
 TGeoArb8.cxx:1676
 TGeoArb8.cxx:1677
 TGeoArb8.cxx:1678
 TGeoArb8.cxx:1679
 TGeoArb8.cxx:1680
 TGeoArb8.cxx:1681
 TGeoArb8.cxx:1682
 TGeoArb8.cxx:1683
 TGeoArb8.cxx:1684
 TGeoArb8.cxx:1685
 TGeoArb8.cxx:1686
 TGeoArb8.cxx:1687
 TGeoArb8.cxx:1688
 TGeoArb8.cxx:1689
 TGeoArb8.cxx:1690
 TGeoArb8.cxx:1691
 TGeoArb8.cxx:1692
 TGeoArb8.cxx:1693
 TGeoArb8.cxx:1694
 TGeoArb8.cxx:1695
 TGeoArb8.cxx:1696
 TGeoArb8.cxx:1697
 TGeoArb8.cxx:1698
 TGeoArb8.cxx:1699
 TGeoArb8.cxx:1700
 TGeoArb8.cxx:1701
 TGeoArb8.cxx:1702
 TGeoArb8.cxx:1703
 TGeoArb8.cxx:1704
 TGeoArb8.cxx:1705
 TGeoArb8.cxx:1706
 TGeoArb8.cxx:1707
 TGeoArb8.cxx:1708
 TGeoArb8.cxx:1709
 TGeoArb8.cxx:1710
 TGeoArb8.cxx:1711
 TGeoArb8.cxx:1712
 TGeoArb8.cxx:1713
 TGeoArb8.cxx:1714
 TGeoArb8.cxx:1715
 TGeoArb8.cxx:1716
 TGeoArb8.cxx:1717
 TGeoArb8.cxx:1718
 TGeoArb8.cxx:1719
 TGeoArb8.cxx:1720
 TGeoArb8.cxx:1721
 TGeoArb8.cxx:1722
 TGeoArb8.cxx:1723
 TGeoArb8.cxx:1724
 TGeoArb8.cxx:1725
 TGeoArb8.cxx:1726
 TGeoArb8.cxx:1727
 TGeoArb8.cxx:1728
 TGeoArb8.cxx:1729
 TGeoArb8.cxx:1730
 TGeoArb8.cxx:1731
 TGeoArb8.cxx:1732
 TGeoArb8.cxx:1733
 TGeoArb8.cxx:1734
 TGeoArb8.cxx:1735
 TGeoArb8.cxx:1736
 TGeoArb8.cxx:1737
 TGeoArb8.cxx:1738
 TGeoArb8.cxx:1739
 TGeoArb8.cxx:1740
 TGeoArb8.cxx:1741
 TGeoArb8.cxx:1742
 TGeoArb8.cxx:1743
 TGeoArb8.cxx:1744
 TGeoArb8.cxx:1745
 TGeoArb8.cxx:1746
 TGeoArb8.cxx:1747
 TGeoArb8.cxx:1748
 TGeoArb8.cxx:1749
 TGeoArb8.cxx:1750
 TGeoArb8.cxx:1751
 TGeoArb8.cxx:1752
 TGeoArb8.cxx:1753
 TGeoArb8.cxx:1754
 TGeoArb8.cxx:1755
 TGeoArb8.cxx:1756
 TGeoArb8.cxx:1757
 TGeoArb8.cxx:1758
 TGeoArb8.cxx:1759
 TGeoArb8.cxx:1760
 TGeoArb8.cxx:1761
 TGeoArb8.cxx:1762
 TGeoArb8.cxx:1763
 TGeoArb8.cxx:1764
 TGeoArb8.cxx:1765
 TGeoArb8.cxx:1766
 TGeoArb8.cxx:1767
 TGeoArb8.cxx:1768
 TGeoArb8.cxx:1769
 TGeoArb8.cxx:1770
 TGeoArb8.cxx:1771
 TGeoArb8.cxx:1772
 TGeoArb8.cxx:1773
 TGeoArb8.cxx:1774
 TGeoArb8.cxx:1775
 TGeoArb8.cxx:1776
 TGeoArb8.cxx:1777
 TGeoArb8.cxx:1778
 TGeoArb8.cxx:1779
 TGeoArb8.cxx:1780
 TGeoArb8.cxx:1781
 TGeoArb8.cxx:1782
 TGeoArb8.cxx:1783
 TGeoArb8.cxx:1784
 TGeoArb8.cxx:1785
 TGeoArb8.cxx:1786
 TGeoArb8.cxx:1787
 TGeoArb8.cxx:1788
 TGeoArb8.cxx:1789
 TGeoArb8.cxx:1790
 TGeoArb8.cxx:1791
 TGeoArb8.cxx:1792
 TGeoArb8.cxx:1793
 TGeoArb8.cxx:1794
 TGeoArb8.cxx:1795
 TGeoArb8.cxx:1796
 TGeoArb8.cxx:1797
 TGeoArb8.cxx:1798
 TGeoArb8.cxx:1799
 TGeoArb8.cxx:1800
 TGeoArb8.cxx:1801
 TGeoArb8.cxx:1802
 TGeoArb8.cxx:1803
 TGeoArb8.cxx:1804
 TGeoArb8.cxx:1805
 TGeoArb8.cxx:1806
 TGeoArb8.cxx:1807
 TGeoArb8.cxx:1808
 TGeoArb8.cxx:1809
 TGeoArb8.cxx:1810
 TGeoArb8.cxx:1811
 TGeoArb8.cxx:1812
 TGeoArb8.cxx:1813
 TGeoArb8.cxx:1814
 TGeoArb8.cxx:1815
 TGeoArb8.cxx:1816
 TGeoArb8.cxx:1817
 TGeoArb8.cxx:1818
 TGeoArb8.cxx:1819
 TGeoArb8.cxx:1820
 TGeoArb8.cxx:1821
 TGeoArb8.cxx:1822
 TGeoArb8.cxx:1823
 TGeoArb8.cxx:1824
 TGeoArb8.cxx:1825
 TGeoArb8.cxx:1826
 TGeoArb8.cxx:1827
 TGeoArb8.cxx:1828
 TGeoArb8.cxx:1829
 TGeoArb8.cxx:1830
 TGeoArb8.cxx:1831
 TGeoArb8.cxx:1832
 TGeoArb8.cxx:1833
 TGeoArb8.cxx:1834
 TGeoArb8.cxx:1835
 TGeoArb8.cxx:1836
 TGeoArb8.cxx:1837
 TGeoArb8.cxx:1838
 TGeoArb8.cxx:1839
 TGeoArb8.cxx:1840
 TGeoArb8.cxx:1841
 TGeoArb8.cxx:1842
 TGeoArb8.cxx:1843
 TGeoArb8.cxx:1844
 TGeoArb8.cxx:1845
 TGeoArb8.cxx:1846
 TGeoArb8.cxx:1847
 TGeoArb8.cxx:1848
 TGeoArb8.cxx:1849
 TGeoArb8.cxx:1850
 TGeoArb8.cxx:1851
 TGeoArb8.cxx:1852
 TGeoArb8.cxx:1853