]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliMagWrapCheb.h
Renamed output file to Vertex.Performance.root
[u/mrichter/AliRoot.git] / STEER / AliMagWrapCheb.h
index 195dc50d71b8899ca24bddf2ff3ce9ce699c964e..6c79770a1f0e51c16b278d0522c07cd30c2951cf 100644 (file)
@@ -5,9 +5,9 @@
 //                                                                               //
 //  Wrapper for the set of mag.field parameterizations by Chebyshev polinomials  //
 //  To obtain the field in cartesian coordinates/components use                  //
-//    Field(float* xyz, float* bxyz);                                            //
+//    Field(double* xyz, double* bxyz);                                          //
 //  For cylindrical coordinates/components:                                      //
-//    FieldCyl(float* rphiz, float* brphiz)                                      //
+//    FieldCyl(double* rphiz, double* brphiz)                                    //
 //                                                                               //
 //  The solenoid part is parameterized in the volume  R<500, -550<Z<550 cm       //
 //                                                                               //
@@ -25,9 +25,9 @@
 //                                                                               //
 //  To obtain the field integral in the TPC region from given point to nearest   //
 //  cathod plane (+- 250 cm) use:                                                //
-//  GetTPCInt(float* xyz, float* bxyz);  for Cartesian frame                     //
+//  GetTPCInt(double* xyz, double* bxyz);  for Cartesian frame                   //
 //  or                                                                           //
-//  GetTPCIntCyl(Float_t *rphiz, Float_t *b); for Cylindrical frame              //
+//  GetTPCIntCyl(Double_t *rphiz, Double_t *b); for Cylindrical frame            //
 //                                                                               //
 //                                                                               //
 //  The units are kiloGauss and cm.                                              //
@@ -57,15 +57,17 @@ class AliMagWrapCheb: public TNamed
   virtual void Clear(const Option_t * = "");
   //
   Int_t      GetNParamsSol()                              const {return fNParamsSol;}
-  Int_t      GetNSegZSol()                                const {return fNSegZSol;}
-  Float_t*     GetSegZSol() const {return fSegZSol;}
+  Int_t      GetNSegZSol()                                const {return fNZSegSol;}
+  Float_t*   GetSegZSol()                                 const {return fSegZSol;}
   //
-  Int_t      GetNParamsTPCInt()                           const {return fNParamsTPCInt;}
-  Int_t      GetNSegZTPCInt()                             const {return fNSegZTPCInt;}
+  Int_t      GetNParamsTPCInt()                           const {return fNParamsTPC;}
+  Int_t      GetNSegZTPCInt()                             const {return fNZSegTPC;}
   //
   Int_t      GetNParamsDip()                              const {return fNParamsDip;}
   Int_t      GetNSegZDip()                                const {return fNZSegDip;}
   //
+  Float_t    GetMaxZ()                                    const {return GetMaxZSol();}
+  Float_t    GetMinZ()                                    const {return fParamsDip ? GetMinZDip() : GetMinZSol();}
   //
   Float_t    GetMinZSol()                                 const {return fMinZSol;}
   Float_t    GetMaxZSol()                                 const {return fMaxZSol;}
@@ -74,54 +76,49 @@ class AliMagWrapCheb: public TNamed
   Float_t    GetMinZDip()                                 const {return fMinZDip;}
   Float_t    GetMaxZDip()                                 const {return fMaxZDip;}
   //
-  Float_t    GetMinZTPCInt()                              const {return fMinZTPCInt;}
-  Float_t    GetMaxZTPCInt()                              const {return fMaxZTPCInt;}
-  Float_t    GetMaxRTPCInt()                              const {return fMaxRTPCInt;}
+  Float_t    GetMinZTPCInt()                              const {return fMinZTPC;}
+  Float_t    GetMaxZTPCInt()                              const {return fMaxZTPC;}
+  Float_t    GetMaxRTPCInt()                              const {return fMaxRTPC;}
   //
   AliCheb3D* GetParamSol(Int_t ipar)                      const {return (AliCheb3D*)fParamsSol->UncheckedAt(ipar);}
-  AliCheb3D* GetParamTPCInt(Int_t ipar)                   const {return (AliCheb3D*)fParamsTPCInt->UncheckedAt(ipar);}
+  AliCheb3D* GetParamTPCInt(Int_t ipar)                   const {return (AliCheb3D*)fParamsTPC->UncheckedAt(ipar);}
   AliCheb3D* GetParamDip(Int_t ipar)                      const {return (AliCheb3D*)fParamsDip->UncheckedAt(ipar);}
   //
   virtual void Print(Option_t * = "")                     const;
   //
-  virtual void Field(const Float_t *xyz, Float_t *b)                const;
-  virtual void Field(const Double_t *xyz, Double_t *b)              const;
+  virtual void Field(const Double_t *xyz, Double_t *b)    const;
+  Double_t     GetBz(const Double_t *xyz)                 const;
   //
-  virtual void FieldCyl(const Float_t *rphiz, Float_t *b)     const;
-  virtual void FieldCyl(const Double_t *rphiz, Double_t *b)   const;
+  void FieldCyl(const Double_t *rphiz, Double_t  *b)      const;
+  void GetTPCInt(const Double_t *xyz, Double_t *b)        const;
+  void GetTPCIntCyl(const Double_t *rphiz, Double_t *b)   const;
   //
-  virtual void GetTPCInt(const Float_t *xyz, Float_t *b)        const;
-  virtual void GetTPCIntCyl(const Float_t *rphiz, Float_t *b)   const;
-  //
-  template <class T>
-    Int_t      FindDipSegment(const T *xyz)               const; 
-  //
-  template <class T>
-    static void CylToCartCylB(const T *rphiz, const T *brphiz,T *bxyz);
-  template <class T>
-    static void CylToCartCartB(const T *xyz,  const T *brphiz,T *bxyz);
-  template <class T>
-    static void CartToCylCartB(const T *xyz,  const T *bxyz,  T *brphiz);
-  template <class T>  
-    static void CartToCylCylB(const T *rphiz, const T *bxyz,  T *brphiz);
-  template <class T>
-    static void CartToCyl(const T *xyz,  T *rphiz);
-  template <class T>
-    static void CylToCart(const T *rphiz,T *xyz);
+  Int_t       FindSolSegment(const Double_t *xyz)         const; 
+  Int_t       FindTPCSegment(const Double_t *xyz)         const; 
+  Int_t       FindDipSegment(const Double_t *xyz)         const; 
+  static void CylToCartCylB(const Double_t *rphiz, const Double_t *brphiz,Double_t *bxyz);
+  static void CylToCartCartB(const Double_t *xyz,  const Double_t *brphiz,Double_t *bxyz);
+  static void CartToCylCartB(const Double_t *xyz,  const Double_t *bxyz,  Double_t *brphiz);
+  static void CartToCylCylB(const Double_t *rphiz, const Double_t *bxyz,  Double_t *brphiz);
+  static void CartToCyl(const Double_t *xyz,  Double_t *rphiz);
+  static void CylToCart(const Double_t *rphiz,Double_t *xyz);
   //
 #ifdef  _INC_CREATION_ALICHEB3D_                          // see AliCheb3D.h for explanation
   void         LoadData(const char* inpfile);
   //
   AliMagWrapCheb(const char* inputFile);
-  void       SaveData(const char* outfile)              const;
-  Int_t      SegmentDipDimension(Float_t** seg,const TObjArray* par,int npar, int dim, 
-                                Float_t xmn,Float_t xmx,Float_t ymn,Float_t ymx,Float_t zmn,Float_t zmx);
+  void       SaveData(const char* outfile)                const;
+  Int_t      SegmentDimension(Float_t** seg,const TObjArray* par,int npar, int dim, 
+                             Float_t xmn,Float_t xmx,Float_t ymn,Float_t ymx,Float_t zmn,Float_t zmx);
   //
   void       AddParamSol(const AliCheb3D* param);
   void       AddParamTPCInt(const AliCheb3D* param);
   void       AddParamDip(const AliCheb3D* param);
-  void       BuildTableDip();
+  void       BuildTable(Int_t npar,TObjArray *parArr, Int_t &nZSeg, Int_t &nYSeg, Int_t &nXSeg,
+                       Float_t &minZ,Float_t &maxZ,Float_t **segZ,Float_t **segY,Float_t **segX,
+                       Int_t **begSegY,Int_t **nSegY,Int_t **begSegX,Int_t **nSegX,Int_t **segID);
   void       BuildTableSol();
+  void       BuildTableDip();
   void       BuildTableTPCInt();
   void       ResetTPCInt();
   //
@@ -129,88 +126,81 @@ class AliMagWrapCheb: public TNamed
 #endif
   //
  protected:
-    virtual void FieldCylSol(const Float_t *rphiz, Float_t *b)      const;
-    virtual void FieldCylSol(const Double_t *rphiz, Double_t *b)    const;
+  void     FieldCylSol(const Double_t *rphiz, Double_t *b)    const;
+  Double_t FieldCylSolBz(const Double_t *rphiz)               const;
   //
  protected:
   //
-  Int_t      fNParamsSol;            // Total number of parameterization pieces for Sol 
-  Int_t      fNSegZSol;              // Number of segments in Z for Solenoid field
-  //
-  Int_t      fNParamsTPCInt;         // Total number of parameterization pieces for TPC field integral 
-  Int_t      fNSegZTPCInt;           // Number of segments in Z for TPC field integral
+  Int_t      fNParamsSol;            // Total number of parameterization pieces for solenoid 
+  Int_t      fNZSegSol;              // number of distinct Z segments in Solenoid
+  Int_t      fNPSegSol;              // number of distinct P segments in Solenoid
+  Int_t      fNRSegSol;              // number of distinct R segments in Solenoid
+  Float_t*   fSegZSol;               //[fNZSegSol] coordinates of distinct Z segments in Solenoid
+  Float_t*   fSegPSol;               //[fNPSegSol] coordinated of P segments for each Zsegment in Solenoid
+  Float_t*   fSegRSol;               //[fNRSegSol] coordinated of R segments for each Psegment in Solenoid
+  Int_t*     fBegSegPSol;            //[fNPSegSol] beginning of P segments array for each Z segment
+  Int_t*     fNSegPSol;              //[fNZSegSol] number of P segments for each Z segment
+  Int_t*     fBegSegRSol;            //[fNPSegSol] beginning of R segments array for each P segment
+  Int_t*     fNSegRSol;              //[fNPSegSol] number of R segments for each P segment
+  Int_t*     fSegIDSol;              //[fNRSegSol] ID of the solenoid parameterization for given RPZ segment
+  Float_t    fMinZSol;               // Min Z of Solenoid parameterization
+  Float_t    fMaxZSol;               // Max Z of Solenoid parameterization
+  TObjArray* fParamsSol;             // Parameterization pieces for Solenoid field
+  Float_t    fMaxRSol;               // max raduis for Solenoid field
+  //
+  Int_t      fNParamsTPC;            // Total number of parameterization pieces for TPCint 
+  Int_t      fNZSegTPC;              // number of distinct Z segments in TPCint
+  Int_t      fNPSegTPC;              // number of distinct P segments in TPCint
+  Int_t      fNRSegTPC;              // number of distinct R segments in TPCint
+  Float_t*   fSegZTPC;               //[fNZSegTPC] coordinates of distinct Z segments in TPCint
+  Float_t*   fSegPTPC;               //[fNPSegTPC] coordinated of P segments for each Zsegment in TPCint
+  Float_t*   fSegRTPC;               //[fNRSegTPC] coordinated of R segments for each Psegment in TPCint
+  Int_t*     fBegSegPTPC;            //[fNPSegTPC] beginning of P segments array for each Z segment
+  Int_t*     fNSegPTPC;              //[fNZSegTPC] number of P segments for each Z segment
+  Int_t*     fBegSegRTPC;            //[fNPSegTPC] beginning of R segments array for each P segment
+  Int_t*     fNSegRTPC;              //[fNPSegTPC] number of R segments for each P segment
+  Int_t*     fSegIDTPC;              //[fNRSegTPC] ID of the TPCint parameterization for given RPZ segment
+  Float_t    fMinZTPC;               // Min Z of TPCint parameterization
+  Float_t    fMaxZTPC;               // Max Z of TPCint parameterization
+  TObjArray* fParamsTPC;             // Parameterization pieces for TPCint field
+  Float_t    fMaxRTPC;               // max raduis for Solenoid field integral in TPC
   //
   Int_t      fNParamsDip;            // Total number of parameterization pieces for dipole 
   Int_t      fNZSegDip;              // number of distinct Z segments in Dipole
   Int_t      fNYSegDip;              // number of distinct Y segments in Dipole
   Int_t      fNXSegDip;              // number of distinct X segments in Dipole
-  //
-  Float_t*   fSegZSol;               //[fNSegZSol]      upper boundaries of Z segments
-  Float_t*   fSegRSol;               //[fNParamsSol]    upper boundaries of R segments
-  //
-  Float_t*   fSegZTPCInt;            //[fNSegZTPCInt]    upper boundaries of Z segments
-  Float_t*   fSegRTPCInt;            //[fNParamsTPCInt]  upper boundaries of R segments
-  //
   Float_t*   fSegZDip;               //[fNZSegDip] coordinates of distinct Z segments in Dipole
   Float_t*   fSegYDip;               //[fNYSegDip] coordinated of Y segments for each Zsegment in Dipole
   Float_t*   fSegXDip;               //[fNXSegDip] coordinated of X segments for each Ysegment in Dipole
-  //
-  Int_t*     fNSegRSol;              //[fNSegZSol]      number of R segments for each Z segment
-  Int_t*     fSegZIdSol;             //[fNSegZSol]      Id of the first R segment of each Z segment in the fSegRSol...
-  //
-  Int_t*     fNSegRTPCInt;           //[fNSegZTPCInt]   number of R segments for each Z segment
-  Int_t*     fSegZIdTPCInt;          //[fNSegZTPCInt]   Id of the first R segment of each Z segment in the fSegRTPCInt...
-  //
   Int_t*     fBegSegYDip;            //[fNZSegDip] beginning of Y segments array for each Z segment
   Int_t*     fNSegYDip;              //[fNZSegDip] number of Y segments for each Z segment
   Int_t*     fBegSegXDip;            //[fNYSegDip] beginning of X segments array for each Y segment
   Int_t*     fNSegXDip;              //[fNYSegDip] number of X segments for each Y segment
   Int_t*     fSegIDDip;              //[fNXSegDip] ID of the dipole parameterization for given XYZ segment
-  //
-  Float_t    fMinZSol;               // Min Z of Sol parameterization (in CYL. coordinates)
-  Float_t    fMaxZSol;               // Max Z of Sol parameterization (in CYL. coordinates)
-  Float_t    fMaxRSol;               // Max R of Sol parameterization (in CYL. coordinates)
-  //
   Float_t    fMinZDip;               // Min Z of Dipole parameterization
   Float_t    fMaxZDip;               // Max Z of Dipole parameterization
-  //
-  Float_t    fMinZTPCInt;            // Min Z of TPCInt parameterization (in CYL. coordinates)
-  Float_t    fMaxZTPCInt;            // Max Z of TPCInt parameterization (in CYL. coordinates)
-  Float_t    fMaxRTPCInt;            // Max R of TPCInt parameterization (in CYL. coordinates)
-  // 
-  TObjArray* fParamsSol;             // Parameterization pieces for Solenoid field
   TObjArray* fParamsDip;             // Parameterization pieces for Dipole field
-  TObjArray* fParamsTPCInt;          // Parameterization pieces for Solenoid field integrals in TPC region
   //
-  ClassDef(AliMagWrapCheb,3)            // Wrapper class for the set of Chebishev parameterizations of Alice mag.field
+  ClassDef(AliMagWrapCheb,5)         // Wrapper class for the set of Chebishev parameterizations of Alice mag.field
   //
  };
 
 
-//__________________________________________________________________________________________
-inline void AliMagWrapCheb::FieldCyl(const Float_t *rphiz, Float_t *b) const
-{
-  // compute field in Cylindircal coordinates
-  //  if (rphiz[2]<GetMinZSol() || rphiz[2]>GetMaxZSol() || rphiz[0]>GetMaxRSol()) {for (int i=3;i--;) b[i]=0; return;}
-  FieldCylSol(rphiz,b);
-}
-
-
 //__________________________________________________________________________________________
 inline void AliMagWrapCheb::FieldCyl(const Double_t *rphiz, Double_t *b) const
 {
   // compute field in Cylindircal coordinates
   //  if (rphiz[2]<GetMinZSol() || rphiz[2]>GetMaxZSol() || rphiz[0]>GetMaxRSol()) {for (int i=3;i--;) b[i]=0; return;}
+  b[0] = b[1] = b[2] = 0;
   FieldCylSol(rphiz,b);
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCartCylB(const T *rphiz, const T *brphiz,T *bxyz)
+inline void AliMagWrapCheb::CylToCartCylB(const Double_t *rphiz, const Double_t *brphiz,Double_t *bxyz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cyl.system
-  T btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
-  T psiPLUSphi = TMath::ATan2(brphiz[1],brphiz[0]) + rphiz[1];
+  Double_t btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
+  Double_t psiPLUSphi = TMath::ATan2(brphiz[1],brphiz[0]) + rphiz[1];
   bxyz[0] = btr*TMath::Cos(psiPLUSphi);
   bxyz[1] = btr*TMath::Sin(psiPLUSphi);
   bxyz[2] = brphiz[2];
@@ -218,12 +208,11 @@ inline void AliMagWrapCheb::CylToCartCylB(const T *rphiz, const T *brphiz,T *bxy
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCartCartB(const T *xyz, const T *brphiz, T *bxyz)
+inline void AliMagWrapCheb::CylToCartCartB(const Double_t* xyz, const Double_t *brphiz, Double_t *bxyz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cart.system
-  T btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
-  T phiPLUSpsi = TMath::ATan2(xyz[1],xyz[0]) +  TMath::ATan2(brphiz[1],brphiz[0]);
+  Double_t btr = TMath::Sqrt(brphiz[0]*brphiz[0]+brphiz[1]*brphiz[1]);
+  Double_t phiPLUSpsi = TMath::ATan2(xyz[1],xyz[0]) +  TMath::ATan2(brphiz[1],brphiz[0]);
   bxyz[0] = btr*TMath::Cos(phiPLUSpsi);
   bxyz[1] = btr*TMath::Sin(phiPLUSpsi);
   bxyz[2] = brphiz[2];
@@ -231,12 +220,11 @@ inline void AliMagWrapCheb::CylToCartCartB(const T *xyz, const T *brphiz, T *bxy
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCylCartB(const T *xyz, const T *bxyz, T *brphiz)
+inline void AliMagWrapCheb::CartToCylCartB(const Double_t *xyz, const Double_t *bxyz, Double_t *brphiz)
 {
   // convert field in cylindrical coordinates to cartesian system, poin is in cart.system
-  T btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
-  T psiMINphi = TMath::ATan2(bxyz[1],bxyz[0]) - TMath::ATan2(xyz[1],xyz[0]);
+  Double_t btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
+  Double_t psiMINphi = TMath::ATan2(bxyz[1],bxyz[0]) - TMath::ATan2(xyz[1],xyz[0]);
   //
   brphiz[0] = btr*TMath::Cos(psiMINphi);
   brphiz[1] = btr*TMath::Sin(psiMINphi);
@@ -245,12 +233,11 @@ inline void AliMagWrapCheb::CartToCylCartB(const T *xyz, const T *bxyz, T *brphi
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCylCylB(const T *rphiz, const T *bxyz, T *brphiz)
+inline void AliMagWrapCheb::CartToCylCylB(const Double_t *rphiz, const Double_t *bxyz, Double_t *brphiz)
 {
   // convert field in cylindrical coordinates to cartesian system, point is in cyl.system
-  T btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
-  T psiMINphi =  TMath::ATan2(bxyz[1],bxyz[0]) - rphiz[1];
+  Double_t btr = TMath::Sqrt(bxyz[0]*bxyz[0]+bxyz[1]*bxyz[1]);
+  Double_t psiMINphi =  TMath::ATan2(bxyz[1],bxyz[0]) - rphiz[1];
   brphiz[0] = btr*TMath::Cos(psiMINphi);
   brphiz[1] = btr*TMath::Sin(psiMINphi);
   brphiz[2] = bxyz[2];
@@ -258,8 +245,7 @@ inline void AliMagWrapCheb::CartToCylCylB(const T *rphiz, const T *bxyz, T *brph
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CartToCyl(const T *xyz,T *rphiz)
+inline void AliMagWrapCheb::CartToCyl(const Double_t *xyz, Double_t *rphiz)
 {
   rphiz[0] = TMath::Sqrt(xyz[0]*xyz[0]+xyz[1]*xyz[1]);
   rphiz[1] = TMath::ATan2(xyz[1],xyz[0]);
@@ -267,32 +253,11 @@ inline void AliMagWrapCheb::CartToCyl(const T *xyz,T *rphiz)
 }
 
 //__________________________________________________________________________________________________
-template <class T>
-inline void AliMagWrapCheb::CylToCart(const T *rphiz, T *xyz)
+inline void AliMagWrapCheb::CylToCart(const Double_t *rphiz, Double_t *xyz)
 {
   xyz[0] = rphiz[0]*TMath::Cos(rphiz[1]);
   xyz[1] = rphiz[0]*TMath::Sin(rphiz[1]);
   xyz[2] = rphiz[2];
 }
 
-//__________________________________________________________________________________________________
-template <class T>
-Int_t    AliMagWrapCheb::FindDipSegment(const T *xyz) const 
-{
-  // find the segment containing point xyz. If it is outside find the closest segment 
-  int xid,yid,zid = TMath::BinarySearch(fNZSegDip,fSegZDip,(Float_t)xyz[2]); // find zsegment
-  int ysegBeg = fBegSegYDip[zid];
-  //
-  for (yid=0;yid<fNSegYDip[zid];yid++) if (xyz[1]<fSegYDip[ysegBeg+yid]) break;
-  if ( --yid < 0 ) yid = 0;
-  yid +=  ysegBeg;
-  //
-  int xsegBeg = fBegSegXDip[yid];
-  for (xid=0;xid<fNSegXDip[yid];xid++) if (xyz[0]<fSegXDip[xsegBeg+xid]) break;
-  if ( --xid < 0) xid = 0;
-  xid +=  xsegBeg;
-  //
-  return fSegIDDip[xid];
-}
-
 #endif