ITS/UPGRADE/KMCDetector.h

   1 #ifndef DETECTORK_H\r
   2 #define DETECTORK_H\r
   3 \r
   4 #include <TNamed.h>\r
   5 #include <TList.h>\r
   6 #include <TGraph.h>\r
   7 #include <Riostream.h>\r
   8 #include <TMatrixD.h>\r
   9 #include <TClonesArray.h>\r
  10 #include <TArrayD.h>\r
  11 #include <TH2.h>\r
  12 #include <TH1.h>\r
  13 #include "AliExternalTrackParam.h"\r
  14 #include "AliLog.h"\r
  15 \r
  16 //-------------------------------------------------------------------------\r
  17 // Current support and development: Ruben Shahoyan (Ruben.Shahoyan@cern.ch) \r
  18 //-------------------------------------------------------------------------\r
  19 \r
  20 \r
  21 class KMCLayer;\r
  22 class KMCCluster;\r
  23 \r
  24 //////////////////////////////////////////////////////////////////////////////////////////////\r
  25 //--------------------------------------------------------------------------------------------\r
  26 class KMCProbe : public AliExternalTrackParam {\r
  27  public:\r
  28   enum {kBitKilled=BIT(14)};\r
  29   enum {kNDOF=5};\r
  30   enum {kY2=0,kZ2=2,kSnp2=5,kTgl2=9,kPtI2=14};\r
  31   enum {kY,kZ,kSnp,kTgl,kPtI};\r
  32   //\r
  33  KMCProbe() : fMass(0.14),fChi2(0),fHits(0),fFakes(0),fNHits(0),fNHitsITS(0),fNHitsITSFake(0),fInnLrCheck(fgNITSLayers) {for (int i=kMaxITSLr;i--;) fClID[i] =-2;\r
  34 }\r
  35   KMCProbe(KMCProbe& src);\r
  36   KMCProbe& operator=(const KMCProbe& src);\r
  37   virtual ~KMCProbe() {}\r
  38   virtual   void   Print(Option_t* option = "") const;\r
  39   virtual   void   Reset() {fMass=0.14; fChi2=0; fHits=fFakes=0;  for (int i=kMaxITSLr;i--;) fClID[i]=-2; AliExternalTrackParam::Reset();}\r
  40   virtual   Bool_t IsSortable()                     const {return kTRUE;}\r
  41   virtual   Int_t  Compare(const TObject* obj)      const;\r
  42   void      ResetCovMat();\r
  43   //\r
  44   void      Kill(Bool_t v=kTRUE)                          {SetBit(kBitKilled,v);}\r
  45   Bool_t    IsKilled()                              const {return TestBit(kBitKilled);}\r
  46   //\r
  47   Bool_t    CorrectForMeanMaterial(const KMCLayer* lr, Bool_t inward=kTRUE);\r
  48   Bool_t    GetXatLabR(Double_t r,Double_t &x, Double_t bz, Int_t dir=0) const;\r
  49   Bool_t    PropagateToR(double r, double b, int dir=0);\r
  50   \r
  51   //\r
  52   void      SetMass(double m=0.14)                      {fMass = m;}\r
  53   Double_t  GetMass()                             const {return fMass;}\r
  54   Double_t  GetChi2()                             const {return fChi2;}\r
  55   void      SetChi2(double chi2)                        {fChi2 = chi2;}\r
  56   void      AddChi2(double chi2)                        {fChi2 += chi2;}\r
  57   void      SetInnerLrChecked(Int_t n)                  {if (n<fgNITSLayers) fInnLrCheck = n;}\r
  58   UShort_t  GetNITSHits()                         const {return fNHitsITS;}\r
  59   UShort_t  GetNFakeITSHits()                     const {return fNHitsITSFake;}\r
  60   UShort_t  GetNHits()                            const {return fNHits;}\r
  61   UShort_t  GetInnerLayerChecked()                const {return fInnLrCheck;}\r
  62   UInt_t&   GetHitsPatt()                               {return fHits;}\r
  63   UInt_t&   GetFakesPatt()                              {return fFakes;}\r
  64   Double_t  GetNormChi2(Bool_t penalize=kFALSE)   const;\r
  65   void      AddHit(Int_t lr, double chi2, Int_t clID=-1);\r
  66   void      ResetHit(Int_t lr);\r
  67   Bool_t    IsHit(Int_t lr)                       const {return (lr<fgNITSLayers) ? IsWBit(fHits,lr)  : kFALSE;}\r
  68   Bool_t    IsHitFake(Int_t lr)                   const {return (lr<fgNITSLayers) ? IsWBit(fFakes,lr) : kFALSE;}\r
  69   Bool_t    PropagateToCluster(KMCCluster* cl, Double_t b);\r
  70   // protected: \r
  71   static void   SetWBit(UInt_t &patt,UInt_t bit)               {patt |= 0x1<<bit;}\r
  72   static void   ResetWBit(UInt_t &patt,UInt_t bit)             {patt &= ~(0x1<<bit);}\r
  73   static Bool_t IsWBit(const UInt_t &patt,const UInt_t bit)    {return patt&(0x1<<bit);}\r
  74   static void   SetNITSLayers(Int_t n)                         {fgNITSLayers = n;}\r
  75   static int    GetNITSLayers()                                {return fgNITSLayers;}\r
  76   //\r
  77   static Double_t GetMissingHitPenalty()                        {return fgMissingHitPenalty;}\r
  78   static void     SetMissingHitPenalty(double p=2.)             {fgMissingHitPenalty = p;}\r
  79   //\r
  80  public:\r
  81   enum {kMaxITSLr=12};\r
  82   Float_t fMass;   // mass\r
  83   Float_t fChi2;   // total chi2\r
  84   UInt_t  fHits;   // pattern on hits (max 32!)\r
  85   UInt_t  fFakes;  // pattern of fakes among hits\r
  86   UShort_t fNHits;    // total hits\r
  87   UShort_t fNHitsITS; // total ITS hits\r
  88   UShort_t fNHitsITSFake; // number of fake ITS hits\r
  89   UShort_t fInnLrCheck;   // lowest active layer where update was checked\r
  90   Int_t    fClID[kMaxITSLr];\r
  91   //\r
  92   static Int_t    fgNITSLayers;\r
  93   static Double_t fgMissingHitPenalty;\r
  94   ClassDef(KMCProbe,2);  \r
  95 };\r
  96 \r
  97 //_______________________________________\r
  98 inline Double_t KMCProbe::GetNormChi2(Bool_t penalize) const\r
  99 {\r
 100   // normalized chi2, penilized for missing hits\r
 101   //  if (fNHitsITS<3) return 0;\r
 102   double chi2 = fChi2;\r
 103   if (penalize) {\r
 104     int nMiss = fgNITSLayers - fNHitsITS - fInnLrCheck;\r
 105     chi2 = fChi2 + fgMissingHitPenalty*nMiss;\r
 106   }\r
 107   return chi2/( (fNHitsITS<<1)-kNDOF);\r
 108 }\r
 109 \r
 110 //_______________________________________\r
 111 inline void KMCProbe::AddHit(Int_t lr, double chi2, Int_t clID) {\r
 112   // note: lr is active layer ID\r
 113   if (lr<0) return;\r
 114   fNHits++;\r
 115   if (lr<fgNITSLayers) {\r
 116     fChi2 += chi2;\r
 117     SetWBit(fHits,lr); \r
 118     fNHitsITS++;\r
 119     if (clID>-1) {\r
 120       SetWBit(fFakes,lr);\r
 121       fNHitsITSFake++;\r
 122     }\r
 123     fClID[lr] = clID;\r
 124     //else ResetWBit(fFakes,lr);\r
 125   }\r
 126 }\r
 127 \r
 128 //_______________________________________\r
 129 inline void KMCProbe::ResetHit(Int_t lr) {\r
 130   // note: lr is active layer ID\r
 131   if (lr>=fgNITSLayers) return; \r
 132   if (IsWBit(fHits,lr))  {fNHitsITS--;     ResetWBit(fHits,lr);}\r
 133   if (IsWBit(fFakes,lr)) {fNHitsITSFake--; ResetWBit(fFakes,lr);}\r
 134 }\r
 135 \r
 136 //////////////////////////////////////////////////////////////////////////////////////////////\r
 137 //--------------------------------------------------------------------------------------------\r
 138 class KMCCluster : public TObject {\r
 139  public:\r
 140   //\r
 141   enum {kBitKilled=BIT(14)};\r
 142  KMCCluster(Float_t y=0, Float_t z=0, Float_t x=0, Float_t phi=0) : fY(y),fZ(z),fX(x),fPhi(phi) {}\r
 143   KMCCluster(KMCCluster &src);\r
 144   KMCCluster& operator=(const KMCCluster& src);\r
 145   virtual ~KMCCluster() {}\r
 146   void Reset()               {Clear();}\r
 147   //\r
 148   Double_t GetY()    const   {return fY;}\r
 149   Double_t GetX()    const   {return fX;}\r
 150   Double_t GetZ()    const   {return fZ;}\r
 151   Double_t GetPhi()  const   {return fPhi;}\r
 152   //\r
 153   void    Kill(Bool_t v=kTRUE)          {SetBit(kBitKilled,v);}\r
 154   Bool_t  IsKilled()              const {return TestBit(kBitKilled);}\r
 155   Float_t fY; \r
 156   Float_t fZ; \r
 157   Float_t fX;\r
 158   Float_t fPhi;\r
 159   void Set(Float_t y, Float_t z, Float_t x, Float_t phi) {fY=y; fZ=z; fX=x; fPhi=phi; ResetBit(kBitKilled);}\r
 160   virtual void Print(Option_t * = 0) const;\r
 161   ClassDef(KMCCluster,1);\r
 162 };\r
 163 \r
 164 //////////////////////////////////////////////////////////////////////////////////////////////\r
 165 //--------------------------------------------------------------------------------------------\r
 166 class KMCLayer : public TNamed {\r
 167 public:\r
 168   enum {kTypeNA=-1,kITS,kTPC,kBitVertex=BIT(15)};\r
 169   KMCLayer(char *name);\r
 170   Float_t GetRadius()   const {return fR;}\r
 171   Float_t GetRadL()     const {return fx2X0;}\r
 172   Float_t GetXTimesRho() const {return fXRho;}\r
 173   Float_t GetPhiRes()   const {return fPhiRes;}\r
 174   Float_t GetZRes()     const {return fZRes;}\r
 175   Float_t GetLayerEff() const {return fEff;}\r
 176   Int_t   GetActiveID() const {return fActiveID;}\r
 177   virtual void  Print(Option_t* option = "") const;\r
 178   //\r
 179   Bool_t IsDead()       const {return fIsDead;}\r
 180   Bool_t IsITS()        const {return fType==kITS;}\r
 181   Bool_t IsTPC()        const {return fType==kTPC;}\r
 182   Bool_t IsVertex()     const {return TestBit(kBitVertex);}\r
 183   //\r
 184   Int_t    AddBgCluster(double y,double z,double x,double phi) {int n=GetNBgClusters(); new (fClBg[n]) KMCCluster(y,z,x,phi); return ++n;}\r
 185   KMCCluster* GetBgCluster(Int_t i) {return (KMCCluster*)fClBg[i];}\r
 186   KMCCluster* GetMCCluster()        {return (KMCCluster*)&fClMC;}\r
 187   //\r
 188   KMCProbe*       GetAnProbe()    const {return (KMCProbe*)&fTrCorr;}\r
 189   Int_t           GetNMCTracks()  const {return fTrMC.GetEntries();}\r
 190   KMCProbe*       AddMCTrack(KMCProbe* src=0);\r
 191   KMCProbe*       GetMCTrack(Int_t it) const {return (KMCProbe*)fTrMC[it];}\r
 192   KMCProbe*       GetWinnerMCTrack()  {if (!fTrMC.IsSorted()) fTrMC.Sort(); return fTrMC.GetEntries() ? (KMCProbe*)fTrMC[0]:0;}\r
 193   TClonesArray*         GetMCTracks()        const {return (TClonesArray*)&fTrMC;}\r
 194   void Reset();\r
 195   void ResetMC() { fClBg.Clear(); fTrMC.Clear();}\r
 196   void ResetBgClusters() { fClBg.Clear(); }\r
 197   void ResetMCTracks()   { fTrMC.Clear(); }\r
 198   Int_t GetNBgClusters() const { return fClBg.GetEntries();}\r
 199   static Double_t GetDefEff()   {return fgDefEff;}\r
 200   static void     SetDefEff(double eff=1) {fgDefEff = eff>1. ? 1.: (eff<0? 0:eff);}\r
 201   //\r
 202   //\r
 203   Float_t fR; \r
 204   Float_t fx2X0;\r
 205   Float_t fXRho;    // x*density\r
 206   Float_t fPhiRes; \r
 207   Float_t fZRes;   \r
 208   Float_t fEff;\r
 209   Bool_t  fIsDead;\r
 210   Int_t   fType;   // its, tpc etc\r
 211   Int_t   fActiveID;   // active layer id\r
 212   Float_t fSig2EstD;\r
 213   Float_t fSig2EstZ;\r
 214   //\r
 215   KMCCluster   fClCorr;     // ideal cluster\r
 216   KMCCluster   fClMC;       // MC cluster (from MS scattered track)\r
 217   TClonesArray fClBg; // bg clusters for MC\r
 218   //\r
 219   KMCProbe     fTrCorr;   // ideal track\r
 220   TClonesArray fTrMC;      // MC tracks\r
 221   //\r
 222   static Double_t fgDefEff;\r
 223   ClassDef(KMCLayer,1);\r
 224 };\r
 225 \r
 226 //////////////////////////////////////////////////////////////////////////////////////////////\r
 227 //--------------------------------------------------------------------------------------------\r
 228 class KMCDetector : public TNamed {\r
 229  public:\r
 230   enum {kUtilHisto=BIT(14)};\r
 231   KMCDetector();\r
 232   KMCDetector(char *name,char *title);\r
 233   virtual ~KMCDetector();\r
 234 \r
 235   void AddLayer(char *name, Float_t radius, Float_t radL, Float_t xrho=0, Float_t phiRes=999999, Float_t zRes=999999, Float_t eff=-1);\r
 236   Int_t GetLayerID(Int_t actID) const;\r
 237   void KillLayer(char *name);\r
 238   void SetRadius(char *name, Float_t radius);\r
 239   void SetRadiationLength(char *name, Float_t radL);\r
 240   void SetResolution(char *name, Float_t phiRes=999999, Float_t zRes=999999);\r
 241   void SetLayerEfficiency(char *name, Float_t eff=1.0);\r
 242   void RemoveLayer(char *name);\r
 243 \r
 244   Float_t GetRadius(char *name);\r
 245   Float_t GetRadiationLength(char *name);\r
 246   Float_t GetResolution(char *name, Int_t axis=0);\r
 247   Float_t GetLayerEfficiency(char *name);\r
 248 \r
 249   void PrintLayout(); \r
 250   void PlotLayout(Int_t plotDead = kTRUE);\r
 251   \r
 252   void MakeAliceAllNew(Bool_t flagTPC =0,  Bool_t flagMon=1,int setVer=0);\r
 253   void MakeAliceCurrent(Bool_t flagTPC =0, Int_t AlignResiduals = 0);\r
 254   void AddTPC(Float_t phiResMean, Float_t zResMean, Int_t skip=1);\r
 255   void RemoveTPC();\r
 256 \r
 257   void SetBField(Float_t bfield) {fBFieldG = bfield*10; }\r
 258   Float_t GetBField() const {return fBFieldG/10; }\r
 259   void SetLhcUPCscale(Float_t lhcUPCscale) {fLhcUPCscale = lhcUPCscale; }\r
 260   Float_t GetLhcUPCscale() const { return fLhcUPCscale; }\r
 261   void SetParticleMass(Float_t particleMass) {fParticleMass = particleMass; }\r
 262   Float_t GetParticleMass() const { return fParticleMass; }\r
 263   void SetIntegrationTime(Float_t integrationTime) {fIntegrationTime = integrationTime; }\r
 264   Float_t GetIntegrationTime() const { return fIntegrationTime; }\r
 265   void SetMaxRadiusOfSlowDetectors(Float_t maxRadiusSlowDet) {fMaxRadiusSlowDet =  maxRadiusSlowDet; }\r
 266   Float_t GetMaxRadiusOfSlowDetectors() const { return fMaxRadiusSlowDet; }\r
 267   void SetAvgRapidity(Float_t avgRapidity) {fAvgRapidity = avgRapidity; }\r
 268   Float_t GetAvgRapidity() const { return fAvgRapidity; }\r
 269   void SetConfidenceLevel(Float_t confLevel) {fConfLevel = confLevel; }\r
 270   Float_t GetConfidenceLevel() const { return fConfLevel; }\r
 271   void SetAtLeastCorr(Int_t atLeastCorr ) {fAtLeastCorr = atLeastCorr; }\r
 272   Float_t GetAtLeastCorr() const { return fAtLeastCorr; }\r
 273   void SetAtLeastFake(Int_t atLeastFake ) {fAtLeastFake = atLeastFake; }\r
 274   Float_t GetAtLeastFake() const { return fAtLeastFake; }\r
 275 \r
 276   void SetdNdEtaCent(Int_t dNdEtaCent ) {fdNdEtaCent = dNdEtaCent; }\r
 277   Float_t GetdNdEtaCent() const { return fdNdEtaCent; }\r
 278 \r
 279   Int_t GetNLayers()          const {return fLayers.GetEntries(); }\r
 280   Int_t GetNActiveLayers()    const {return fNActiveLayers; }\r
 281   Int_t GetNActiveITSLayers() const {return fNActiveITSLayers; }\r
 282 \r
 283   Bool_t SolveSingleTrackViaKalman(Double_t mass, Double_t pt, Double_t eta);\r
 284   Bool_t SolveSingleTrackViaKalmanMC(int offset=6);\r
 285   Bool_t SolveSingleTrack(Double_t mass, Double_t pt, Double_t eta, TObjArray* sumArr=0, int nMC=10000,int offset=6);\r
 286   KMCProbe* KalmanSmooth(int actLr, int actMin,int actMax) const;\r
 287   KMCProbe* KalmanSmoothFull(int actLr, int actMin,int actMax) const; //TBD\r
 288   void   EliminateUnrelated();\r
 289   //\r
 290   KMCProbe* PrepareKalmanTrack(double pt, double lambda, double mass, int charge, double phi=0,double x=0,double y=0,double z=0);\r
 291   Bool_t TransportKalmanTrackWithMS(KMCProbe *probTr, int maxLr) const;\r
 292   void   ApplyMS(KMCProbe* trc,  double x2x0) const;\r
 293   Bool_t PropagateToLayer(KMCProbe* trc, KMCLayer* lr, int dir) const;\r
 294   Bool_t UpdateTrack(KMCProbe* trc, KMCLayer* lr, KMCCluster* cl, Bool_t goToCluster=kTRUE) const;\r
 295 \r
 296   // Helper functions\r
 297   Double_t ThetaMCS                 ( Double_t mass, Double_t RadLength, Double_t momentum ) const;\r
 298   Double_t ProbGoodHit              ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; \r
 299   Double_t ProbGoodChiSqHit         ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; \r
 300   Double_t ProbGoodChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; \r
 301   Double_t ProbNullChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; \r
 302  \r
 303   // Howard W. hit distribution and convolution integral\r
 304   Double_t Dist              ( Double_t Z, Double_t radius ) ;  \r
 305   Double_t HitDensity        ( Double_t radius )   ;\r
 306   Double_t UpcHitDensity     ( Double_t radius )   ;\r
 307   Double_t IntegratedHitDensity  ( Double_t multiplicity, Double_t radius )   ;\r
 308   Double_t OneEventHitDensity    ( Double_t multiplicity, Double_t radius ) const   ;\r
 309   Double_t D0IntegratedEfficiency( Double_t pt, Double_t corrEfficiency[][20] ) const ;\r
 310   \r
 311   TGraph* GetGraphMomentumResolution(Int_t color, Int_t linewidth=1);\r
 312   TGraph* GetGraphPointingResolution(Int_t axis,Int_t color, Int_t linewidth=1);\r
 313   TGraph* GetGraphPointingResolutionTeleEqu(Int_t axis,Int_t color, Int_t linewidth=1);\r
 314 \r
 315   TGraph* GetGraphImpactParam(Int_t mode, Int_t axis, Int_t color, Int_t linewidth=1);\r
 316 \r
 317   TGraph* GetGraphRecoEfficiency(Int_t particle, Int_t color, Int_t linewidth=1); \r
 318   TGraph* GetGraphRecoFakes(Int_t particle,Int_t color, Int_t linewidth);\r
 319 \r
 320   TGraph* GetGraph(Int_t number, Int_t color, Int_t linewidth=1);\r
 321 \r
 322   void MakeStandardPlots(Bool_t add =0, Int_t color=1, Int_t linewidth=1,Bool_t onlyPionEff=0);\r
 323 \r
 324   // method to extend AliExternalTrackParam functionality\r
 325   Bool_t IsZero(double val, double tol=1e-9) const {return TMath::Abs(val)<tol;}\r
 326   TList *GetLayers()                   const {return (TList*)&fLayers;}\r
 327   KMCLayer* GetLayer(Int_t i)          const {return (KMCLayer*) fLayers.At(i);}\r
 328   KMCLayer* GetActiveLayer(Int_t actID)    const {int pid=GetLayerID(actID); return pid<0 ? 0:GetLayer(pid);}\r
 329   KMCLayer* GetLayer(const char* name) const {return (KMCLayer*) fLayers.FindObject(name);}\r
 330   KMCProbe* GetProbeTrack()       const {return (KMCProbe*)&fProbe;}\r
 331   void   ClassifyLayers();\r
 332   void   ResetMCTracks(int maxLr=-1);                      \r
 333   //\r
 334   Bool_t   GetUseBackground()               const {return fUseBackground;}\r
 335   void     SetUseBackground(Bool_t v=kTRUE)       {fUseBackground = v;}\r
 336   void     CheckTrackProlongations(KMCProbe *probe, KMCLayer* lr, KMCLayer* lrP);\r
 337   void     ResetSearchLimits() {fBgYMin=fBgZMin=1e6; fBgYMax=fBgZMax=-1e6; fNBgLimits=0;}\r
 338   void     UpdateSearchLimits(KMCProbe* probe, KMCLayer* lr);\r
 339   Int_t    GenBgClusters(KMCLayer* lr);\r
 340   Bool_t   NeedToKill(KMCProbe* probe) const;\r
 341   Double_t PropagateBack(KMCProbe* trc);\r
 342   //\r
 343   // definition of reconstructable track\r
 344   void     RequireMaxChi2Cl(double cut=25.)           {fMaxChi2Cl = cut>0 ? cut:9; fMaxChi2ClSQ = TMath::Sqrt(fMaxChi2Cl);}\r
 345   void     RequireMinITSHits(Int_t n=4)               {fMinITSHits = n;}\r
 346   void     RequireMaxNormChi2NDF(double cut=5.)       {fMaxNormChi2NDF = cut>0 ? cut:9;}\r
 347   void     RequirePattern(UInt_t *patt, int groups);\r
 348   //\r
 349   Double_t GetMaxChi2Cl()                      const {return fMaxChi2Cl;}\r
 350   Double_t GetMaxNormChi2NDFusterKMC()              const {return fMaxNormChi2NDF;}\r
 351   Int_t    GetMinITSHits()                     const {return fMinITSHits;}\r
 352   //\r
 353   Double_t GetUpdCalls()                       const {return fUpdCalls;}\r
 354   void     InitMCWatchHistos();\r
 355   TH2F*    GetHMCLrResidRPhi()                 const {return fHMCLrResidRPhi;}\r
 356   TH2F*    GetHMCLrResidZ()                    const {return fHMCLrResidZ;}\r
 357   TH2F*    GetHMCLrChi2()                      const {return fHMCLrChi2;}\r
 358   //\r
 359   void     PrintITS(Option_t* opt="") const {for (int i=0;i<=fLastActiveITSLayer;i++) if (!GetLayer(i)->IsDead()) GetLayer(i)->Print(opt);}\r
 360   static void SetVtxConstraint(double d=-1, double z=-1) {fgVtxConstraint[0]=d; fgVtxConstraint[1]=z;}\r
 361   //\r
 362   void CalcHardSearchLimits(double dzv);\r
 363   void SetMaxSeedToPropagate(Int_t n=300) {fMaxSeedToPropagate = n;}\r
 364 \r
 365  protected:\r
 366  \r
 367   Int_t fNLayers;        // total number of layers in the model\r
 368   Int_t fNActiveLayers;  // number of active layers in the model\r
 369   Int_t fNActiveITSLayers;  // number of active ITS layers in the model\r
 370   Int_t fLastActiveLayer;       // id of last active layer\r
 371   Int_t fLastActiveITSLayer;    // id of last active ITS layer\r
 372   Int_t fLastActiveLayerTracked;    // id of last active layer really used for tracking of given pt\r
 373   TList fLayers;                // List of layer pointers\r
 374   Float_t fBFieldG;             // Magnetic Field in Gauss (set in Tesla)\r
 375   Float_t fLhcUPCscale;         // UltraPeripheralElectrons: scale from RHIC to LHC\r
 376   Float_t fIntegrationTime;     // electronics integration time\r
 377   Float_t fConfLevel;           // Confidence Level for the tracking\r
 378   Float_t fAvgRapidity;         // rapidity of the track (= mean)\r
 379   Float_t fParticleMass;        // Particle used for tracking. Standard: mass of pion\r
 380   Double_t fMaxRadiusSlowDet;   // Maximum radius for slow detectors.  Fast detectors \r
 381                                 // and only fast detectors reside outside this radius.\r
 382   Int_t fAtLeastCorr;     // min. number of correct hits for the track to be "good"\r
 383   Int_t fAtLeastFake;     // min. number of fake hits for the track to be "fake"\r
 384 \r
 385   Int_t fdNdEtaCent;       // Multiplicity\r
 386   //\r
 387   // reconstruction settings\r
 388   Double_t fMaxChi2Cl;   // max cluster-track chi2 \r
 389   Double_t fMaxNormChi2NDF;// max chi2/NDF to accept\r
 390   Int_t    fMinITSHits;  // min ITS hits in track to accept\r
 391   //\r
 392   Double_t fMaxChi2ClSQ; // precalculated sqrt(chi2);\r
 393   //\r
 394   Int_t    fMaxSeedToPropagate; // take 1st fMaxSeedToPropagate seeds from each layer\r
 395   // background settings\r
 396   Bool_t   fUseBackground; // do we want to simulate background?\r
 397   Double_t fBgYMin;      // min Y for current layer bg generation\r
 398   Double_t fBgYMax;      // max Y for current layer bg generation\r
 399   Double_t fBgZMin;      // min Z for current layer bg generation\r
 400   Double_t fBgZMax;      // max Z for current layer bg generation\r
 401   TArrayD  fBgYMinTr;    // min Y for each seed at current layer\r
 402   TArrayD  fBgYMaxTr;    // max Y for each seed at current layer\r
 403   TArrayD  fBgZMinTr;    // min Z for each seed at current layer\r
 404   TArrayD  fBgZMaxTr;    // max Z for each seed at current layer\r
 405   Int_t    fNBgLimits;   // depth of limits array\r
 406   //\r
 407   enum {kMaxNDetectors = 200};\r
 408  \r
 409   Double_t fTransMomenta[40];                          // array of transverse momenta\r
 410   Double_t fMomentumRes[40];                           // array of momentum resolution\r
 411   Double_t fResolutionRPhi[40];                        // array of rphi resolution\r
 412   Double_t fResolutionZ[40];                           // array of z resolution\r
 413   Double_t fDetPointRes[kMaxNDetectors][40];           // array of rphi resolution per layer\r
 414   Double_t fDetPointZRes[kMaxNDetectors][40];          // array of z resolution per layer\r
 415   Double_t fEfficiency[3][40];                         // efficiency for different particles\r
 416   Double_t fFake[3][40];                               // fake prob for different particles\r
 417   //\r
 418   KMCProbe fProbe;\r
 419   Double_t fDensFactorEta;                             // density scaling for non-0 eta (precalculated\r
 420   //\r
 421   Double_t fUpdCalls;                  // number of kalman updates\r
 422   TH2F*  fHMCLrResidRPhi;              // Residials on lr, rphi\r
 423   TH2F*  fHMCLrResidZ;                 // Residials on lr, rphi\r
 424   TH2F*  fHMCLrChi2;                   // track to clusyer chi2 on each lr, rphi\r
 425   TArrayI fPattITS;                     // bit pattern of each group of active layers where hit is required\r
 426   //\r
 427   static Double_t fgVtxConstraint[2];  // if both positive, the vertex is used as constraint (accounted in chi2 but not in update)\r
 428   ClassDef(KMCDetector,1);\r
 429 };\r
 430 \r
 431 //____________________________________________________________________________\r
 432 inline Bool_t KMCProbe::PropagateToCluster(KMCCluster* cl, double b)\r
 433 {\r
 434   // propagate track to cluster frame\r
 435   if (!Rotate(cl->GetPhi()) || !PropagateTo(cl->GetX(),b)) {\r
 436     AliDebug(2,Form("Failed to propager track to cluster at phi=%.3f X=%.3f",cl->GetPhi(),cl->GetX()));\r
 437     if (AliLog::GetGlobalDebugLevel()>1) Print();\r
 438     return kFALSE;\r
 439   }\r
 440   return kTRUE;\r
 441 }\r
 442 \r
 443 //////////////////////////////////////////////////////////////////////////////////////////////\r
 444 //--------------------------------------------------------------------------------------------\r
 445 //\r
 446 // Class to collect summary of certrain track performance\r
 447 // The track can be selected according to min and max ITS clusters to accept (total, fake and or correct)\r
 448 // and according to presence of (any and/or correct)  or absence of fakes clusters in certain groups\r
 449 // of ITS layers.\r
 450 \r
 451 class KMCTrackSummary : public TNamed\r
 452 {\r
 453  public:\r
 454   KMCTrackSummary(const char* name=0,const char* title=0,int nlr=0);\r
 455   virtual ~KMCTrackSummary();\r
 456   //\r
 457   void Add(const KMCTrackSummary* src, double scl=1);\r
 458 \r
 459   virtual void  Print(Option_t* option = "") const;\r
 460   KMCTrackSummary* MakeCopy(const char* pref) const;\r
 461   void   SetNamePrefix(const char* pref="");\r
 462 \r
 463   Bool_t CheckTrack(KMCProbe* trc);\r
 464   void   AddTrack(KMCProbe* trc);\r
 465   //\r
 466   TH1*   GetHMCChi2()               const {return fHMCChi2;}\r
 467   TH1*   GetHMCSigDCARPhi()         const {return fHMCSigDCARPhi;}\r
 468   TH1*   GetHMCSigDCAZ()            const {return fHMCSigDCAZ;}\r
 469   TH1*   GetHMCSigPt()              const {return fHMCSigPt;}\r
 470   TH2*   GetHMCNCl()                const {return fHMCNCl;}\r
 471   TH1*   GetHMCFakePatt()           const {return fHMCFakePatt;}\r
 472   //\r
 473   Double_t GetCountAcc()            const {return fCountAcc;}\r
 474   Double_t GetCountTot()            const {return fCountTot;}\r
 475   Double_t GetEff()                 const {return fCountTot>0 ? fCountAcc/fCountTot : 0.;}\r
 476   Double_t GetEffErr()              const;\r
 477   KMCProbe& GetAnProbe()      const {return (KMCProbe&)fAnProbe;}\r
 478   KMCProbe& GetRefProbe()     const {return (KMCProbe&)fRefProbe;}\r
 479   void SetAnProbe(KMCProbe* pr)     {if (pr) fAnProbe = *pr;}\r
 480   void SetRefProbe(KMCProbe* pr)    {if (pr) fRefProbe = *pr;}\r
 481   //\r
 482   void SetMinMaxClITS(Int_t mn=0,Int_t mx=999)     {fMinMaxClITS[0]=mn; fMinMaxClITS[1]=mx;}\r
 483   void SetMinMaxClITSFake(Int_t mn=0,Int_t mx=999) {fMinMaxClITSFake[0]=mn; fMinMaxClITSFake[1]=mx;}\r
 484   void SetMinMaxClITSCorr(Int_t mn=0,Int_t mx=999) {fMinMaxClITSCorr[0]=mn; fMinMaxClITSCorr[1]=mx;}\r
 485   //\r
 486   Double_t GetUpdCalls()           const {return fCountTot>0 ? fUpdCalls/fCountTot : 0;}\r
 487   void     AddUpdCalls(double v)         {fUpdCalls += v;}\r
 488   //\r
 489   void AddPatternITS(Int_t patt)           {AddPattern(patt,fPattITS);}\r
 490   void AddPatternITSFakeExcl(Int_t patt)   {AddPattern(patt,fPattITSFakeExcl);}\r
 491   void AddPatternITSCorr(Int_t patt)       {AddPattern(patt,fPattITSCorr);}\r
 492   static UInt_t Bits(Bool_t l0=0, Bool_t l1=0, Bool_t l2=0, Bool_t l3=0, Bool_t l4=0, Bool_t l5=0, Bool_t l6=0, Bool_t l7=0,Bool_t  l8=0, Bool_t l9=0,\r
 493                      Bool_t l10=0,Bool_t l11=0,Bool_t l12=0,Bool_t l13=0,Bool_t l14=0,Bool_t l15=0,Bool_t l16=0,Bool_t l17=0,Bool_t l18=0,Bool_t l19=0,\r
 494                      Bool_t l20=0,Bool_t l21=0,Bool_t l22=0,Bool_t l23=0,Bool_t l24=0,Bool_t l25=0,Bool_t l26=0,Bool_t l27=0,Bool_t l28=0,Bool_t l29=0,\r
 495                      Bool_t l30=0,Bool_t l31=0);\r
 496   //\r
 497   // protected:\r
 498   void   AddPattern(Int_t patt,TArrayI& arr) {int n=arr.GetSize(); arr.Set(n+1); arr[n]=patt;}\r
 499   Bool_t CheckPattern(UInt_t patt, Int_t cond) {return ((UInt_t)cond)&patt;}\r
 500   Int_t  OutOfRange(Int_t n, Int_t range[2])  {if (n<range[0]) return -1; if (n>range[1]) return 1; return 0;}\r
 501   static void   PrependTNamed(TNamed* obj, const char* nm=0, const char* tit=0);\r
 502   //\r
 503  private:\r
 504  KMCTrackSummary(const KMCTrackSummary& src) :TNamed(src) {}\r
 505  protected:\r
 506   //\r
 507   Int_t   fNITSLayers;                  // number of its layers\r
 508   // track selection conditions\r
 509   Int_t   fMinMaxClITS[2];              // min and max N ITS clusters (total)\r
 510   Int_t   fMinMaxClITSFake[2];          // min and max N ITS fake clusters  \r
 511   Int_t   fMinMaxClITSCorr[2];          // min and max N ITS correct clusters\r
 512   //\r
 513   TArrayI fPattITS;                     // bit pattern of each group of layers where hit is required\r
 514   TArrayI fPattITSFakeExcl;             // bit pattern of each group of layers where fake hit is not tolerated\r
 515   TArrayI fPattITSCorr;                 // bit pattern of each group of layers where correc hit is required\r
 516   //\r
 517   TH1F*  fHMCChi2;                     // track chi2\r
 518   TH1F*  fHMCSigDCARPhi;               // DCA distribution in RPhi from MC\r
 519   TH1F*  fHMCSigDCAZ;                  // DCA distribution in Z for from MC\r
 520   TH1F*  fHMCSigPt;                    // Pt difference distribution from MC\r
 521   TH2F*  fHMCNCl;                      // number of clusters\r
 522   TH1F*  fHMCFakePatt;                 // Fakes pattern\r
 523   //\r
 524   Double_t fCountAcc;                  // track counter of accepted tracks\r
 525   Double_t fCountTot;                  // track counter of all tracks\r
 526   Double_t fUpdCalls;                  // total update calls\r
 527   //\r
 528   KMCProbe fRefProbe;                  // reference probe\r
 529   KMCProbe fAnProbe;                   // analital probe\r
 530   static Int_t fgSumCounter;           // global counter of bins\r
 531 \r
 532   ClassDef(KMCTrackSummary,1)\r
 533 };\r
 534 \r
 535 inline Double_t KMCTrackSummary::GetEffErr() const {\r
 536   //\r
 537   if (fCountTot<1) return 0;\r
 538   double r = fCountAcc/fCountTot;\r
 539   double err = r*(1.-r)/fCountTot;\r
 540   return err>0 ? TMath::Sqrt(err) : 0;\r
 541 }\r
 542 \r
 543 inline KMCTrackSummary* KMCTrackSummary::MakeCopy(const char* pref) const {\r
 544   // create a copy, prepending all names with prefix\r
 545   KMCTrackSummary* sm = (KMCTrackSummary*) Clone(this->GetName()); \r
 546   sm->SetNamePrefix(pref); \r
 547   return sm;\r
 548 } \r
 549 \r
 550 #endif\r