RICH/AliRICHParam.h

   1 #ifndef AliRICHParam_h
   2 #define AliRICHParam_h
   3
   4 #include <TError.h>
   5 #include <TMath.h>
   6 #include <TObjArray.h>
   7 #include <TObject.h>
   8 #include <TMath.h>
   9 #include <TRandom.h>
  10 #include <TVector.h>
  11 #include <TVector2.h>
  12 #include <TVector3.h>
  13 #include <TRandom.h>
  14 #include <TError.h>
  15 #include <TObjArray.h>
  16 #include <AliLog.h>
  17 #include <TClass.h>
  18
  19 static const int kNchambers=7;     //number of RICH chambers
  20 static const int kNpadsX = 160;    //number of pads along X in single chamber
  21 static const int kNpadsY = 144;    //number of pads along Y in single chamber
  22 static const int kNsectors=6;      //number of sectors per chamber
  23
  24 static const int kCerenkov=50000050;  //??? go to something more general like TPDGCode
  25 static const int kFeedback=50000051;  //??? go to something more general like TPDGCode
  26
  27 class AliRICHChamber;
  28
  29 // Class providing all the needed parametrised information
  30 // to construct the geometry, to define segmentation and to provide response model
  31 // In future will also provide all the staff needed for alignment and calibration
  32
  33
  34 class AliRICHParam :public TObject
  35 {
  36 public:
  37 //ctor&dtor
  38                   AliRICHParam():TObject(),fpChambers(0)  {CreateChambers();}
  39   virtual        ~AliRICHParam()                          {delete fpChambers;}
  40 //test methodes
  41          void     Print(Option_t *opt="") const;                                         //print current parametrization
  42          void     Test()                            {TestSeg();TestTrans();TestResp();}  //test all groups of methodes
  43          void     TestResp();                                                            //test the response group of methodes
  44          void     TestSeg();                                                             //test the segmentation group of methodes
  45          void     TestTrans();                                                           //test the transform group of methodes
  46   static void     DrawAxis();
  47   static void     DrawSectors();
  48 //flags staff
  49   static           void     SetAerogel(Bool_t a)                   {fgIsAerogel=a;}
  50   static           Bool_t   IsAerogel()                            {return fgIsAerogel;}
  51   static           void     SetRadioSrc(Bool_t a)                   {fgIsRadioSrc=a;}
  52   static           Bool_t   IsRadioSrc()                            {return fgIsRadioSrc;}
  53   static              void     SetTestBeam(Bool_t a)                {fgIsTestBeam=a;}
  54   static              Bool_t   IsTestBeam()                         {return fgIsTestBeam;}
  55   static                void     SetWireSag(Bool_t a)               {fgIsWireSag=a;}
  56   static                Bool_t   IsWireSag()                        {return fgIsWireSag;}
  57   static                   void     SetResolveClusters(Bool_t a)    {fgIsResolveClusters=a;}
  58   static                   Bool_t   IsResolveClusters()             {return fgIsResolveClusters;}
  59 //Chambers manipulation  methodes
  60   void            CreateChambers();                                                      //form chamber structure
  61   AliRICHChamber* C(Int_t i)                 {return (AliRICHChamber*)fpChambers->UncheckedAt(i-1);}      //returns pointer to chamber i
  62   Int_t           Nchambers()                {return fpChambers->GetEntriesFast();}      //returns number of chambers
  63 //Geometrical properties
  64   static        Int_t      NpadsX()                   {return kNpadsX;}                           //pads along X in chamber
  65   static        Int_t      NpadsY()                   {return kNpadsY;}                           //pads along Y in chamber
  66   static        Int_t      NpadsXsec()                {return NpadsX()/2;}                        //pads along X in sector
  67   static        Int_t      NpadsYsec()                {return NpadsY()/3;}                        //pads along Y in sector
  68   static        Double_t   DeadZone()                 {return 2.6;}                               //dead zone size in cm
  69   static        Double_t   SectorSizeX()              {return NpadsX()*PadSizeX()/2;}             //sector size x, cm
  70   static        Double_t   SectorSizeY()              {return NpadsY()*PadSizeY()/3;}             //sector size y, cm
  71   static        Double_t   PcSizeX()                  {return NpadsX()*PadSizeX()+DeadZone();}    //PC size x, cm
  72   static        Double_t   PcSizeY()                  {return NpadsY()*PadSizeY()+2*DeadZone();}  //PC size y, cm
  73   static        Double_t   Zfreon()                   {return 1.5;}                               //freon thinkness, cm
  74   static        Double_t   Zwin()                     {return 0.5;}                               //radiator quartz window, cm
  75   static        Double_t   Pc2Win()                   {return 8.0;}                               //cm between CsI PC and radiator quartz window
  76   static        Double_t   Pc2Coll()                  {return 7.0;}                               //cm between CsI PC and third wire grid (collection wires)
  77   static        Double_t   Pc2Anod()                  {return 0.204;}                             //cm between CsI PC and first wire grid (anod wires)
  78   static        Double_t   Pc2Cath()                  {return 0.445;}                             //cm between CsI PC and second wire grid (cathode wires)
  79   static        Double_t   Freon2Pc()                 {return Zfreon()+Zwin()+Pc2Win();}          //cm between CsI PC and entrance to freon
  80   static        Double_t   PitchAnod()                {return PadSizeY()/2;}                      //cm between anode wires
  81   static        Double_t   PitchCath()                {return PadSizeY()/4;}                      //cm between cathode wires
  82   static        Double_t   PitchColl()                {return 0.5;}                               //cm between collection wires
  83   static        Double_t   PadSizeX()                 {return 0.8;}                               //pad size x, cm
  84   static        Double_t   PadSizeY      (                               ){return 0.84;}                              //pad size y, cm
  85 //trasformation methodes
  86   static        Int_t      Pad2Cha       (Int_t pad                      ){return pad/100000000;                     }//abs pad -> chamber
  87   static        Int_t      Pad2Sec       (Int_t pad                      ){return pad%100000000/1000000;             }//abs pad -> sector
  88   static        Int_t      Pad2PadX      (Int_t pad                      ){return pad%1000000/1000;                  }//abs pad -> pad x
  89   static        Int_t      Pad2PadY      (Int_t pad                      ){return pad%1000000%100;                   }//abs pad -> pad y
  90   static        Int_t      PadAbs        (Int_t c,Int_t s,Int_t x,Int_t y){return 100000000*c+1000000*s+1000*x+y;    }//(c,s,x,y) -> abs pad
  91   static inline TVector2   Pad2Loc       (Int_t pad                      );                                           //abs pad ->LORS
  92   static inline TVector2   Pad2Loc       (TVector pad                    );                                           //pad  -> LORS returns center of the pad
  93   static        TVector2   Pad2Loc       (Int_t x,Int_t y                ){TVector pad(2);pad[0]=x;pad[1]=y;return Pad2Loc(pad);}//return center of the pad (x,y)
  94   static inline TVector    Loc2Area      (const TVector2 &x2             );                                           //pads area affected by hit x2. area is LeftDown-RightUp pad numbers
  95   static inline Int_t      Loc2Sec       (const TVector2 &x2             );                                           //LORS -> sector
  96   static        Int_t      Loc2Sec       (Double_t x,Double_t y          ){return Loc2Sec(TVector2(x,y));}            //LORS -> sector
  97   static inline TVector    Loc2Pad       (const TVector2 &x2             );                                           //LORS -> pad
  98   static        TVector    Loc2Pad       (Double_t x,Double_t y          ){return Loc2Pad(TVector2(x,y));}            //LORS -> pad
  99   static inline Int_t      Pad2Sec       (const TVector &pad             );                                           //pad  -> sector
 100   static inline Int_t      PadNeighbours (Int_t iPadX,Int_t iPadY,Int_t aListX[4],Int_t aListY[4]);                   //pad -> list of it neighbours
 101   static        Bool_t     IsAccepted    (const TVector2 &x2             ){return ( x2.X()>=0 && x2.X()<=PcSizeX() && x2.Y()>=0 && x2.Y()<=PcSizeY() ) ? kTRUE:kFALSE;}
 102 //optical properties methodes
 103   static        Double_t   MeanCkovEnergy(                               ){return 6.766;}                 //mean Ckov energy according to the total trasmission curve
 104   static        Float_t    PhotonEnergy  (Int_t i                        ){return 0.1*i+5.5;}                         //photon energy (eV) for i-th point
 105   static        Float_t    AbsCH4        (Float_t ev                     );                                           //CH4 abs len (cm)
 106   static        Float_t    AbsGel        (Float_t                        ){return 500;}                               //Aerogel abs len (cm)
 107   static        Float_t    RefIdxC6F14   (Float_t eV                     ){return eV*0.0172+1.177;}                   //Freon ref idx
 108   static        Float_t    RefIdxCH4     (Float_t                        ){return 1.000444;}                          //Methane ref idx
 109   static        Float_t    RefIdxSiO2    (Float_t eV                     ){Float_t e1=10.666,e2=18.125,f1=46.411,f2= 228.71; return TMath::Sqrt(1.+f1/(e1*e1-eV*eV)+f2/(e2*e2-eV*eV));}//Quartz window ref index from TDR p.35
 110   static        Float_t    RefIdxGel     (Float_t                        ){return 1.05;}                              //aerogel ref index
 111   static        Float_t    DenGel        (                               ){return (RefIdxGel(0)-1)/0.21;}             //aerogel density gr/cm^3 parametrization by E.Nappi
 112
 113
 114   static Double_t IonisationPotential()      {return 26.0e-9;}                            //for CH4 in GeV taken from ????
 115   static TVector2 MathiesonDelta()           {return TVector2(5*0.18,5*0.18);}            //area of 5 sigmas of Mathieson distribution (cm)
 116   static Int_t    MaxQdc()                   {return 4095;}                               //QDC number of channels
 117
 118   static Int_t    QthMIP()                   {return 100;}
 119   static Double_t DmatchMIP()                {return 1.;}
 120   static Double_t PmodMax()                  {return 6.5;}
 121   static Int_t    HV(Int_t sector)           {if (sector>=1 && sector <=6) return fgHV[sector-1];  else return -1;} //high voltage for this sector
 122   static void     SetHV(Int_t sector,Int_t hv){fgHV[sector-1]=hv;}
 123 //charge response methodes
 124   inline static Double_t Mathieson(Double_t x1,Double_t x2,Double_t y1,Double_t y2);               //Mathienson integral over given limits
 125   inline static Double_t GainSag(Double_t x,Int_t sector);                                         //gain variations in %
 126          static Double_t QdcSlope(Int_t sec){switch(sec){case -1: return 0;  default:   return 33;}} //weight of electon in QDC channels
 127          static Double_t Gain(const TVector2 &x2){//gives chamber gain in terms of QDC channels for given point in local ref system
 128                           if(fgIsWireSag) return QdcSlope(Loc2Sec(x2))*(1+GainSag(x2.X(),Loc2Sec(x2))/100);
 129                           else            return QdcSlope(Loc2Sec(x2));}
 130   inline static Double_t FracQdc(const TVector2 &x2,const TVector &pad);                           //charge fraction to pad from hit
 131   inline static Int_t    TotQdc(TVector2 x2,Double_t eloss);                                       //total charge for Eloss (GeV) 0 for photons
 132   inline static Bool_t   IsOverTh(Int_t c,TVector pad,Double_t q);                                 //is QDC of the pad registered by FEE
 133          static Int_t    NsigmaTh()                    {return fgNsigmaTh;}                        //
 134          static Float_t  SigmaThMean()                 {return fgSigmaThMean;}                     //QDC electronic noise mean
 135          static Float_t  SigmaThSpread()               {return fgSigmaThSpread;}                   //QDC electronic noise width
 136
 137          static Double_t CogCorr(Double_t x) {return 3.31267e-2*TMath::Sin(2*TMath::Pi()/PadSizeX()*x) //correction of cluster CoG due to sinoidal
 138                                                     -2.66575e-3*TMath::Sin(4*TMath::Pi()/PadSizeX()*x)
 139                                                     +2.80553e-3*TMath::Sin(6*TMath::Pi()/PadSizeX()*x)+0.0070;}
 140          static void     ReadErrFiles();                                                                  //Read Err file parameters
 141          static TVector3 SigmaSinglePhoton(Int_t Npart, Double_t mom, Double_t theta, Double_t phi);      //Find Sigma for single photon from momentum and particle id
 142          static TVector3 SigmaSinglePhoton(Double_t thetaCer, Double_t theta, Double_t phi);              //Fing sigma for single photon from thetacer
 143          static Double_t Interpolate(Double_t par[4][330],Double_t x, Double_t y, Double_t phi);          //Find the error value from interpolation
 144
 145          static TVector3 ForwardTracing(TVector3 entranceTrackPoint,TVector3 vectorTrack, Double_t thetaC, Double_t phiC); //it traces foward a photon from Emission Point to PC
 146          static TVector3 PlaneIntersect(TVector3 vstart,TVector3 p0,TVector3 n,TVector3 v0);              //it finds intersection between straight track and plane
 147          static Double_t SnellAngle(Float_t n1, Float_t n2, Float_t theta1);                              // Snell law
 148          static void     AnglesInDRS(Double_t trackTheta,Double_t trackPhi,Double_t thetaCerenkov,Double_t phiCerenkov,Double_t &tout,Double_t &pout);//It finds photon angles in
 149                                                                                                                                                       //Detector Reference System
 150
 151   static Bool_t     fgIsAerogel;                            //aerogel geometry instead of normal RICH flag
 152   static Double_t fgMass[5];                                // mass array
 153 protected:
 154   static Bool_t     fgIsRadioSrc;                           //radioactive source instead of radiators flag
 155   static Bool_t     fgIsTestBeam;                           //test beam geometry instead of normal RICH flag
 156   static Bool_t     fgIsWireSag;                            //wire sagitta ON/OFF flag
 157   static Bool_t     fgIsResolveClusters;                    //declustering ON/OFF flag
 158   static Bool_t     fgIsFeedback;                           //generate feedback photon?
 159
 160          TObjArray *fpChambers;                             //list of chambers
 161   static Int_t      fgHV[6];                                //HV applied to anod wires
 162   static Int_t      fgNsigmaTh;                             //n. of sigmas to cut for zero suppression
 163   static Float_t    fgSigmaThMean;                          //sigma threshold value
 164   static Float_t    fgSigmaThSpread;                        //spread of sigma
 165
 166   static Double_t fgErrChrom[4][330];                       //
 167   static Double_t fgErrGeom[4][330];                        //
 168   static Double_t fgErrLoc[4][330];                         //Chromatic, Geometric and Localization array to parametrize SigmaCerenkov
 169
 170   ClassDef(AliRICHParam,6)                                  //RICH main parameters class
 171 };
 172 //__________________________________________________________________________________________________
 173 Int_t AliRICHParam::PadNeighbours(Int_t iPadX,Int_t iPadY,Int_t listX[4],Int_t listY[4])
 174 {
 175 //Determines all the neighbouring pads for the given one (iPadX,iPadY). Returns total number of these pads.
 176 //Dead zones are taken into account, meaning pads from different sector are not taken.
 177 //   1
 178 // 2   3
 179 //   4
 180   Int_t nPads=0;
 181   if(iPadY!=NpadsY()&&iPadY!=2*NpadsYsec()&&iPadY!=NpadsYsec()){listX[nPads]=iPadX;   listY[nPads]=iPadY+1; nPads++;}       //1
 182   if(iPadX!=1&&iPadX!=NpadsXsec()+1)                           {listX[nPads]=iPadX-1; listY[nPads]=iPadY;   nPads++;}       //2
 183   if(iPadX!=NpadsXsec()&&iPadX!=NpadsX())                      {listX[nPads]=iPadX+1; listY[nPads]=iPadY;   nPads++;}       //3
 184   if(iPadY!=1&&iPadY!=NpadsYsec()+1&&2*NpadsYsec()+1)          {listX[nPads]=iPadX;   listY[nPads]=iPadY-1; nPads++;}       //4
 185
 186   return nPads;
 187 }//Pad2ClosePads()
 188 //__________________________________________________________________________________________________
 189 Int_t AliRICHParam::Loc2Sec(const TVector2 &v2)
 190 {
 191 //Determines sector containing the given point.
 192 //Returns sector code:
 193 //y ^  5 6
 194 //  |  3 4
 195 //  |  1 2
 196 //   -------> x
 197   Double_t x0=0; Double_t x1=SectorSizeX(); Double_t x2=SectorSizeX()+DeadZone(); Double_t x3=PcSizeX();
 198   Double_t y0=0; Double_t y1=SectorSizeY(); Double_t y2=SectorSizeY()+DeadZone(); Double_t y3=2*SectorSizeY()+DeadZone();
 199   Double_t y4=PcSizeY()-SectorSizeY();      Double_t y5=PcSizeY();
 200
 201   Int_t sector=-1;
 202   if     (v2.X() >= x0 && v2.X() <= x1 )  sector=1;
 203   else if(v2.X() >= x2 && v2.X() <= x3 )  sector=2;
 204   else                                    return -1;
 205
 206   if     (v2.Y() >= y0 && v2.Y() <= y1 )  ;                    //sectors 1 or 2
 207   else if(v2.Y() >= y2 && v2.Y() <= y3 )  sector+=2;           //sectors 3 or 4
 208   else if(v2.Y() >= y4 && v2.Y() <= y5 )  sector+=4;           //sectors 5 or 6
 209   else                                    return -1;
 210   return sector;
 211 }//Loc2Sec(Double_t x, Double_t y)
 212 //__________________________________________________________________________________________________
 213 TVector AliRICHParam::Loc2Pad(const TVector2 &loc)
 214 {
 215 //Determines pad number TVector(padx,pady) containing the given point x2 defined in the chamber RS.
 216 //Pad count starts in lower left corner from 1,1 to 144,160 in upper right corner of a chamber.
 217 //y ^  5 6
 218 //  |  3 4
 219 //  |  1 2
 220 //   -------> x
 221   TVector pad(2);
 222   Int_t sec=Loc2Sec(loc);//trasforms x2 to sector reference system
 223   if(sec==-1) {pad[0]=pad[1]=-1; return pad;}
 224 //first we deal with x
 225   if(sec==1||sec==3||sec==5)    pad[0]=           Int_t(            loc.X()   / PadSizeX() )+1; //sector 1 or 3 or 5
 226   else                          pad[0]=NpadsX() - Int_t( (PcSizeX()-loc.X())  / PadSizeX() )  ; //sector 2 or 4 or 6
 227 //second deal with y
 228        if(sec==1||sec==2)       pad[1]=Int_t(             loc.Y()                / PadSizeY())+1;               //sector 1 or 2
 229   else if(sec==3||sec==4)       pad[1]=Int_t( (loc.Y()-SectorSizeY()-DeadZone()) / PadSizeY())+NpadsYsec()+1;  //sector 3 or 4
 230   else                          pad[1]=NpadsY() - Int_t( (PcSizeY()-loc.Y())     / PadSizeY());                //sector 5 or 6
 231   return pad;
 232 }
 233 //__________________________________________________________________________________________________
 234 Int_t AliRICHParam::Pad2Sec(const TVector &pad)
 235 {
 236 //Determines sector containing the given pad.
 237   Int_t sector=-1;
 238   if     (pad[0] >= 1           && pad[0] <=   NpadsXsec() )    {sector=1;}
 239   else if(pad[0] >  NpadsXsec() && pad[0] <=   NpadsX()    )    {sector=2;}
 240   else                                                         AliDebugClass(1,Form("Wrong pad (%3.0f,%3.0f)",pad[0],pad[1]));
 241
 242   if     (pad[1] >= 1             && pad[1] <=   NpadsYsec() )    {}
 243   else if(pad[1] >  NpadsYsec()   && pad[1] <= 2*NpadsYsec() )    {sector+=2;}
 244   else if(pad[1] >  2*NpadsYsec() && pad[1] <=   NpadsY()    )    {sector+=4;}
 245   else                                                         AliDebugClass(1,Form("Wrong pad (%3.0f,%3.0f)",pad[0],pad[1]));
 246
 247   return sector;
 248 }//Pad2Sec()
 249 //__________________________________________________________________________________________________
 250 TVector2 AliRICHParam::Pad2Loc(TVector pad)
 251 {
 252 //Returns position of the center of the given pad in local system of the chamber (cm)
 253 // y ^  5 6
 254 //   |  3 4        sector numbers
 255 //   |  1 2
 256 //    -------> x
 257   Double_t x=-1,y=-1;
 258   if(pad[0] > 0 && pad[0] <= NpadsXsec())//it's 1 or 3 or 5
 259     x=(pad[0]-0.5)*PadSizeX();
 260   else if(pad[0] > NpadsXsec() && pad[0] <= NpadsX())//it's 2 or 4 or 6
 261     x=(pad[0]-0.5)*PadSizeX()+DeadZone();
 262   else
 263     AliDebugClass(1,Form("Wrong pad (%3.0f,%3.0f)",pad[0],pad[1]));
 264
 265   if(pad[1] > 0 && pad[1] <= NpadsYsec())//it's 1 or 2
 266     y=(pad[1]-0.5)*PadSizeY();
 267   else if(pad[1] > NpadsYsec() && pad[1] <= 2*NpadsYsec())//it's 3 or 4
 268     y=(pad[1]-0.5)*PadSizeY()+DeadZone();
 269   else if(pad[1] > 2*NpadsYsec() && pad[1]<= NpadsY())//it's 5 or 6
 270     y=(pad[1]-0.5)*PadSizeY()+2*DeadZone();
 271   else
 272     AliDebugClass(1,Form("Wrong pad (%3.0f,%3.0f)",pad[0],pad[1]));
 273
 274   return TVector2(x,y);
 275 }
 276 //__________________________________________________________________________________________________
 277 TVector2 AliRICHParam::Pad2Loc(Int_t pad)
 278 {
 279 //Converts absolute pad number to local position in LORS
 280 //LORS is a chamber  reference system with origin in left-down coner looking from IP
 281 //Arguments: pad- absolute pad number
 282 //  Returns: pad center position as TVector2 in PCRS
 283   TVector2 pos;
 284   pos.Set((Pad2PadX(pad)-0.5)*PadSizeX() , (Pad2PadY(pad)-0.5)*PadSizeY());//set to sector LORS
 285   return pos;
 286 }
 287 //__________________________________________________________________________________________________
 288 Double_t AliRICHParam::GainSag(Double_t x,Int_t sector)
 289 {
 290 //Returns % of gain variation due to wire sagita.
 291 //All curves are parametrized as per sector basis, so x must be apriory transformed to the Sector RS.
 292 //Here x is a distance along wires.
 293   x-=SectorSizeX()/2;
 294   if(x>SectorSizeX()) x-=SectorSizeX();
 295   switch(HV(sector)){
 296     case 2150: return 9e-6*TMath::Power(x,4)+2e-7*TMath::Power(x,3)-0.0316*TMath::Power(x,2)-3e-4*x+25.367;//%
 297     case 2100: return 8e-6*TMath::Power(x,4)+2e-7*TMath::Power(x,3)-0.0283*TMath::Power(x,2)-2e-4*x+23.015;
 298     case 2050: return 7e-6*TMath::Power(x,4)+1e-7*TMath::Power(x,3)-0.0254*TMath::Power(x,2)-2e-4*x+20.888;
 299     case 2000: return 6e-6*TMath::Power(x,4)+8e-8*TMath::Power(x,3)-0.0227*TMath::Power(x,2)-1e-4*x+18.961;
 300     default:   return 0;
 301   }
 302 }
 303 //__________________________________________________________________________________________________
 304 Int_t AliRICHParam::TotQdc(TVector2 x2,Double_t eloss)
 305 {
 306 //Calculates the total charge produced by the eloss in point x2 (Chamber RS).
 307 //Returns this change parametrised in QDC channels, or 0 if the hit in the dead zone.
 308 //eloss=0 means photon which produces 1 electron only eloss > 0 for Mip
 309   if(Loc2Sec(x2)==-1) return 0; //hit in the dead zone
 310   Int_t iNelectrons=Int_t(eloss/IonisationPotential()); if(iNelectrons==0) iNelectrons=1;
 311   Double_t qdc=0;
 312   for(Int_t i=1;i<=iNelectrons;i++) qdc+=-Gain(x2)*TMath::Log(gRandom->Rndm());
 313   return Int_t(qdc);
 314 }
 315 //__________________________________________________________________________________________________
 316 Double_t AliRICHParam::FracQdc(const TVector2 &x2,const TVector &pad)
 317 {
 318 //Calculates the charge fraction induced to given pad by the hit from the given point.
 319 //Integrated Mathieson distribution is used.
 320   TVector2 center2=Pad2Loc(pad);//gives center of requested pad
 321   Double_t normXmin=(x2.X()-center2.X()-PadSizeX()/2)  /Pc2Cath();//parametrise for Mathienson
 322   Double_t normXmax=(x2.X()-center2.X()+PadSizeX()/2)  /Pc2Cath();
 323   Double_t normYmin=(x2.Y()-center2.Y()-PadSizeY()/2)  /Pc2Cath();
 324   Double_t normYmax=(x2.Y()-center2.Y()+PadSizeY()/2)  /Pc2Cath();
 325
 326 //requested pad might not belong to the sector of the given hit position, hence the check:
 327   return (Loc2Sec(x2)!=Pad2Sec(pad)) ? 0:Mathieson(normXmin, normYmin, normXmax, normYmax);
 328 }
 329 //__________________________________________________________________________________________________
 330 Double_t AliRICHParam::Mathieson(Double_t x1,Double_t y1,Double_t x2,Double_t y2)
 331 {
 332 //This is the answer to electrostatic problem of charge distrubution in MWPC described elsewhere. (NIM A370(1988)602-603)
 333 //Arguments: x1- diff between center of distribution and left margin of interested pad divided by anod-cathode distance
 334 //           x2,y1,y2- analogically
 335 //  Returns: a charge fraction [0-1].
 336   const Double_t kSqrtKx3=0.77459667;const Double_t kX2=0.962;const Double_t kX4=0.379;
 337   const Double_t kSqrtKy3=0.77459667;const Double_t kY2=0.962;const Double_t kY4=0.379;
 338
 339   Double_t ux1=kSqrtKx3*TMath::TanH(kX2*x1);
 340   Double_t ux2=kSqrtKx3*TMath::TanH(kX2*x2);
 341   Double_t uy1=kSqrtKy3*TMath::TanH(kY2*y1);
 342   Double_t uy2=kSqrtKy3*TMath::TanH(kY2*y2);
 343   return 4*kX4*(TMath::ATan(ux2)-TMath::ATan(ux1))*kY4*(TMath::ATan(uy2)-TMath::ATan(uy1));
 344 }
 345 //__________________________________________________________________________________________________
 346 TVector AliRICHParam::Loc2Area(const TVector2 &x2)
 347 {
 348 //Calculates the area of disintegration for a given point. It's assumed here that this points lays on anode wire.
 349 //Area is a rectangulare set of pads defined by its left-down and right-up coners.
 350   TVector area(4);
 351   TVector pad=Loc2Pad(x2);
 352   area[0]=area[2]=pad[0]; area[1]=area[3]=pad[1];//area is just a pad fired
 353   if(pad[0]!=1           && pad[0]!= NpadsXsec()+1                            ) area[0]--; //left down coner X
 354   if(pad[1]!=1           && pad[1]!= NpadsYsec()+1 && pad[1]!= 2*NpadsYsec()+1) area[1]--; //left down coner Y
 355   if(pad[0]!=NpadsXsec() && pad[0]!= NpadsX()                                 ) area[2]++; //right up coner X
 356   if(pad[1]!=NpadsYsec() && pad[1]!= 2*NpadsYsec() && pad[1]!= NpadsY()       ) area[3]++; //right up coner Y
 357   return area;
 358 }
 359 //__________________________________________________________________________________________________
 360 Bool_t AliRICHParam::IsOverTh(Int_t ,TVector ,Double_t q)
 361 {
 362 //Checks if the current q is over threshold and FEE will save this value to data concentrator.
 363   return (q>NsigmaTh()*(SigmaThMean()+(1.-2*gRandom->Rndm())*SigmaThSpread()));
 364 }
 365 //__________________________________________________________________________________________________
 366 #endif //AliRICHParam_h