[u/mrichter/AliRoot.git] / ITS / UPGRADE / DetectorK.h

#ifndef DETECTORK_H
#define DETECTORK_H

#include <TNamed.h>
#include <TList.h>
#include <TGraph.h>
#include <Riostream.h>

/***********************************************************

Fast Simulation tool for Inner Tracker Systems

original code of using the billoir technique was developed
for the HFT (STAR), James H. Thomas, jhthomas@lbl.gov
http://rnc.lbl.gov/~jhthomas

Changes by S. Rossegger

July 2011  - Adding the possibility of "fake calculation" and "efficiency calculation" 
             with a number of "at least correct clusters on the track"
             Done using the complete combinatorics table with 3^nLayer track outcomes.

April 2011 - Now uses the Kalman method (aliroot implementation) instead of the Billoir
             technique ... (changes by Ruben Shahoyan)

March 2011 - Changes to comply with the Alice Offline coding conventions

Feb. 2011 - Improvement in "lowest pt allowed" -> now uses helix param. for calc. of a,b

          - Adding a more sophisticaed efficiency calculation which includes
            the possibility to make chi2 cuts via Confidence Levels (method of Ruben Shahoyan)
            plus adding 'basic' detection efficiencies per layer ...

          - Adding "ITS Stand alone" tracking capabilities via 
            forward+backward tracking -> Kalman smoothing is then 
            used for the parameter estimates (important for efficiencies)

Jan. 2011 - Inclusion of ImpactParameter Plots (with vtx estimates)
            to allow comparison with ITS performances in pp data

Dec. 2010 - Translation into C++ class format 
          - Adding various Setters and Getters to build the geometry 
	    (based on cylinders) plus handling of the layer properties 


***********************************************************/

class AliExternalTrackParam; 
#include <TMatrixD.h>

class DetectorK : public TNamed {

 public:
  
  DetectorK();
  DetectorK(char *name,char *title);
  virtual ~DetectorK();

  enum {kNptBins = 100}; // less then 400 !!
 
  void AddLayer(char *name, Float_t radius, Float_t radL, Float_t phiRes=999999, Float_t zRes=999999, Float_t eff=0.95);

  void KillLayer(char *name);
  void SetRadius(char *name, Float_t radius);
  void SetRadiationLength(char *name, Float_t radL);
  void SetResolution(char *name, Float_t phiRes=999999, Float_t zRes=999999);
  void SetLayerEfficiency(char *name, Float_t eff=0.95);
  void RemoveLayer(char *name);

  Float_t GetRadius(char *name);
  Float_t GetRadiationLength(char *name);
  Float_t GetResolution(char *name, Int_t axis=0);
  Float_t GetLayerEfficiency(char *name);

  void PrintLayout(); 
  void PlotLayout(Int_t plotDead = kTRUE);
  
  void MakeAliceAllNew(Bool_t flagTPC =1,Bool_t flagMon=1);
  void MakeAliceCurrent(Int_t AlignResiduals = 0, Bool_t flagTPC =1);
  void AddTPC(Float_t phiResMean, Float_t zResMean, Int_t skip=1);
  void RemoveTPC();

  void SetBField(Float_t bfield) {fBField = bfield; }
  Float_t GetBField() const {return fBField; }
  void SetLhcUPCscale(Float_t lhcUPCscale) {fLhcUPCscale = lhcUPCscale; }
  Float_t GetLhcUPCscale() const { return fLhcUPCscale; }
  void SetParticleMass(Float_t particleMass) {fParticleMass = particleMass; }
  Float_t GetParticleMass() const { return fParticleMass; }
  void SetIntegrationTime(Float_t integrationTime) {fIntegrationTime = integrationTime; }
  Float_t GetIntegrationTime() const { return fIntegrationTime; }
  void SetMaxRadiusOfSlowDetectors(Float_t maxRadiusSlowDet) {fMaxRadiusSlowDet =  maxRadiusSlowDet; }
  Float_t GetMaxRadiusOfSlowDetectors() const { return fMaxRadiusSlowDet; }
  void SetAvgRapidity(Float_t avgRapidity) {fAvgRapidity = avgRapidity; }
  Float_t GetAvgRapidity() const { return fAvgRapidity; }
  void SetConfidenceLevel(Float_t confLevel) {fConfLevel = confLevel; }
  Float_t GetConfidenceLevel() const { return fConfLevel; }
  void SetAtLeastCorr(Int_t atLeastCorr ) {fAtLeastCorr = atLeastCorr; }
  Int_t GetAtLeastCorr() const { return fAtLeastCorr; }
  void SetAtLeastFake(Int_t atLeastFake ) {fAtLeastFake = atLeastFake; }
  Int_t GetAtLeastFake() const { return fAtLeastFake; }

  void SetMaxSeedRadius(Double_t maxSeedRadius ) {fMaxSeedRadius = maxSeedRadius; }
  Int_t GetMaxSeedRadius() const { return fMaxSeedRadius; }

  void SetptScale(Double_t ptScale ) {fptScale = ptScale; }
  Int_t GetptScale() const { return fptScale; }


  void SetdNdEtaCent(Int_t dNdEtaCent ) {fdNdEtaCent = dNdEtaCent; }
  Float_t GetdNdEtaCent() const { return fdNdEtaCent; }
  
  
  Float_t GetNumberOfActiveLayers() const {return fNumberOfActiveLayers; }
  Float_t GetNumberOfActiveITSLayers() const {return fNumberOfActiveITSLayers; }

  void SolveViaBilloir(Int_t flagD0=1,Int_t print=1, Bool_t allPt=1, Double_t meanPt =0.250, char* detLayer=((char*)""));

  // Helper functions
  Double_t ThetaMCS                 ( Double_t mass, Double_t RadLength, Double_t momentum ) const;
  Double_t ProbGoodHit              ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; 
  Double_t ProbGoodChiSqHit         ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; 
  Double_t ProbGoodChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; 
  Double_t ProbNullChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ )   ; 
 
  // Howard W. hit distribution and convolution integral
  Double_t Dist              ( Double_t Z, Double_t radius ) ;  
  Double_t HitDensity        ( Double_t radius )   ;
  Double_t UpcHitDensity     ( Double_t radius )   ;
  Double_t IntegratedHitDensity  ( Double_t multiplicity, Double_t radius )   ;
  Double_t OneEventHitDensity    ( Double_t multiplicity, Double_t radius ) const   ;
  Double_t D0IntegratedEfficiency( Double_t pt, Double_t corrEfficiency[][400] ) const ;
  
  TGraph* GetGraphMomentumResolution(Int_t color, Int_t linewidth=1);
  TGraph* GetGraphPointingResolution(Int_t axis,Int_t color, Int_t linewidth=1);
  TGraph* GetGraphPointingResolutionTeleEqu(Int_t axis,Int_t color, Int_t linewidth=1);
  TGraph* GetGraphLayerInfo(Int_t plot, Int_t color, Int_t linewidth=1);

  TGraph* GetGraphImpactParam(Int_t mode, Int_t axis, Int_t color, Int_t linewidth=1);

  TGraph* GetGraphRecoEfficiency(Int_t particle, Int_t color, Int_t linewidth=1); 
  TGraph* GetGraphRecoFakes(Int_t particle,Int_t color, Int_t linewidth);
  TGraph* GetGraphRecoPurity(Int_t particle,Int_t color, Int_t linewidth);

  TGraph* GetGraph(Int_t number, Int_t color, Int_t linewidth=1);

  void MakeStandardPlots(Bool_t add =0, Int_t color=1, Int_t linewidth=1,Bool_t onlyPionEff=0);

  // method to extend AliExternalTrackParam functionality
  static Bool_t GetXatLabR(AliExternalTrackParam* tr,Double_t r,Double_t &x, Double_t bz, Int_t dir=0);
  static Bool_t PropagateToR(AliExternalTrackParam* trc, double r, double b, int dir=0);
  Double_t* PrepareEffFakeKombinations(TMatrixD *probKomb, TMatrixD *probLay);

 protected:
 
  Int_t fNumberOfLayers;        // total number of layers in the model
  Int_t fNumberOfActiveLayers;  // number of active layers in the model
  Int_t fNumberOfActiveITSLayers;  // number of active ITS layers in the model
  TList fLayers;                // List of layer pointers
  Float_t fBField;              // Magnetic Field in Tesla
  Float_t fLhcUPCscale;         // UltraPeripheralElectrons: scale from RHIC to LHC
  Float_t fIntegrationTime;     // electronics integration time
  Float_t fConfLevel;           // Confidence Level for the tracking
  Float_t fAvgRapidity;         // rapidity of the track (= mean)
  Float_t fParticleMass;        // Particle used for tracking. Standard: mass of pion
  Double_t fMaxRadiusSlowDet;   // Maximum radius for slow detectors.  Fast detectors 
                                // and only fast detectors reside outside this radius.
  Int_t fAtLeastCorr;     // min. number of correct hits for the track to be "good"
  Int_t fAtLeastFake;     // min. number of fake hits for the track to be "fake"

  Double_t fMaxSeedRadius; // max seeding radius, e.g. to exclude TRD
  Double_t fptScale;      // used for maxPt

  Int_t fdNdEtaCent;       // Multiplicity

  enum {kMaxNumberOfDetectors = 200};
 
  Double_t fTransMomenta[kNptBins];                          // array of transverse momenta
  Double_t fMomentumRes[kNptBins];                           // array of momentum resolution
  Double_t fResolutionRPhi[kNptBins];                        // array of rphi resolution
  Double_t fResolutionZ[kNptBins];                           // array of z resolution
  Double_t fDetPointRes[kMaxNumberOfDetectors][kNptBins];    // array of rphi resolution per layer
  Double_t fDetPointZRes[kMaxNumberOfDetectors][kNptBins];   // array of z resolution per layer
  Double_t fEfficiency[3][kNptBins];                         // efficiency for different particles
  Double_t fFake[3][kNptBins];                               // fake prob for different particles

  Int_t kDetLayer;                              // layer for which a few more details are extracted
  Double_t fResolutionRPhiLay[kNptBins];                        // array of rphi resolution
  Double_t fResolutionZLay[kNptBins];                           // array of z resolution
  Double_t fEfficProlongLay[kNptBins];                           // array of z resolution

  ClassDef(DetectorK,1);
};

#endif
Commit	Line	Data
84356182	1	#ifndef DETECTORK_H
	2	#define DETECTORK_H
	3
	4	#include <TNamed.h>
	5	#include <TList.h>
	6	#include <TGraph.h>
	7	#include <Riostream.h>
	8
	9	/***********************************************************
	10
	11	Fast Simulation tool for Inner Tracker Systems
	12
	13	original code of using the billoir technique was developed
	14	for the HFT (STAR), James H. Thomas, jhthomas@lbl.gov
	15	http://rnc.lbl.gov/~jhthomas
	16
	17	Changes by S. Rossegger
	18
	19	July 2011 - Adding the possibility of "fake calculation" and "efficiency calculation"
	20	with a number of "at least correct clusters on the track"
	21	Done using the complete combinatorics table with 3^nLayer track outcomes.
	22
	23	April 2011 - Now uses the Kalman method (aliroot implementation) instead of the Billoir
	24	technique ... (changes by Ruben Shahoyan)
	25
	26	March 2011 - Changes to comply with the Alice Offline coding conventions
	27
	28	Feb. 2011 - Improvement in "lowest pt allowed" -> now uses helix param. for calc. of a,b
	29
	30	- Adding a more sophisticaed efficiency calculation which includes
	31	the possibility to make chi2 cuts via Confidence Levels (method of Ruben Shahoyan)
	32	plus adding 'basic' detection efficiencies per layer ...
	33
	34	- Adding "ITS Stand alone" tracking capabilities via
	35	forward+backward tracking -> Kalman smoothing is then
	36	used for the parameter estimates (important for efficiencies)
	37
	38	Jan. 2011 - Inclusion of ImpactParameter Plots (with vtx estimates)
	39	to allow comparison with ITS performances in pp data
	40
	41	Dec. 2010 - Translation into C++ class format
	42	- Adding various Setters and Getters to build the geometry
	43	(based on cylinders) plus handling of the layer properties
	44
	45
	46	***********************************************************/
	47
	48	class AliExternalTrackParam;
	49	#include <TMatrixD.h>
	50
	51	class DetectorK : public TNamed {
	52
	53	public:
	54
	55	DetectorK();
	56	DetectorK(char name,char title);
	57	virtual ~DetectorK();
	58
	59	enum {kNptBins = 100}; // less then 400 !!
	60
	61	void AddLayer(char *name, Float_t radius, Float_t radL, Float_t phiRes=999999, Float_t zRes=999999, Float_t eff=0.95);
	62
	63	void KillLayer(char *name);
	64	void SetRadius(char *name, Float_t radius);
65	void SetRadiationLength(char *name, Float_t radL);
66	void SetResolution(char *name, Float_t phiRes=999999, Float_t zRes=999999);
67	void SetLayerEfficiency(char *name, Float_t eff=0.95);
68	void RemoveLayer(char *name);
69
70	Float_t GetRadius(char *name);
71	Float_t GetRadiationLength(char *name);
72	Float_t GetResolution(char *name, Int_t axis=0);
73	Float_t GetLayerEfficiency(char *name);
74
75	void PrintLayout();
76	void PlotLayout(Int_t plotDead = kTRUE);
77
78	void MakeAliceAllNew(Bool_t flagTPC =1,Bool_t flagMon=1);
79	void MakeAliceCurrent(Int_t AlignResiduals = 0, Bool_t flagTPC =1);
80	void AddTPC(Float_t phiResMean, Float_t zResMean, Int_t skip=1);
81	void RemoveTPC();
82
83	void SetBField(Float_t bfield) {fBField = bfield; }
84	Float_t GetBField() const {return fBField; }
85	void SetLhcUPCscale(Float_t lhcUPCscale) {fLhcUPCscale = lhcUPCscale; }
86	Float_t GetLhcUPCscale() const { return fLhcUPCscale; }
87	void SetParticleMass(Float_t particleMass) {fParticleMass = particleMass; }
88	Float_t GetParticleMass() const { return fParticleMass; }
89	void SetIntegrationTime(Float_t integrationTime) {fIntegrationTime = integrationTime; }
90	Float_t GetIntegrationTime() const { return fIntegrationTime; }
91	void SetMaxRadiusOfSlowDetectors(Float_t maxRadiusSlowDet) {fMaxRadiusSlowDet = maxRadiusSlowDet; }
92	Float_t GetMaxRadiusOfSlowDetectors() const { return fMaxRadiusSlowDet; }
93	void SetAvgRapidity(Float_t avgRapidity) {fAvgRapidity = avgRapidity; }
94	Float_t GetAvgRapidity() const { return fAvgRapidity; }
95	void SetConfidenceLevel(Float_t confLevel) {fConfLevel = confLevel; }
96	Float_t GetConfidenceLevel() const { return fConfLevel; }
97	void SetAtLeastCorr(Int_t atLeastCorr ) {fAtLeastCorr = atLeastCorr; }
98	Int_t GetAtLeastCorr() const { return fAtLeastCorr; }
99	void SetAtLeastFake(Int_t atLeastFake ) {fAtLeastFake = atLeastFake; }
100	Int_t GetAtLeastFake() const { return fAtLeastFake; }
101
102	void SetMaxSeedRadius(Double_t maxSeedRadius ) {fMaxSeedRadius = maxSeedRadius; }
103	Int_t GetMaxSeedRadius() const { return fMaxSeedRadius; }
104
105	void SetptScale(Double_t ptScale ) {fptScale = ptScale; }
106	Int_t GetptScale() const { return fptScale; }
107
108
109
110	void SetdNdEtaCent(Int_t dNdEtaCent ) {fdNdEtaCent = dNdEtaCent; }
111	Float_t GetdNdEtaCent() const { return fdNdEtaCent; }
112
113
114
115	Float_t GetNumberOfActiveLayers() const {return fNumberOfActiveLayers; }
116	Float_t GetNumberOfActiveITSLayers() const {return fNumberOfActiveITSLayers; }
117
118	void SolveViaBilloir(Int_t flagD0=1,Int_t print=1, Bool_t allPt=1, Double_t meanPt =0.250, char* detLayer=((char*)""));
119
120	// Helper functions
121	Double_t ThetaMCS ( Double_t mass, Double_t RadLength, Double_t momentum ) const;
122	Double_t ProbGoodHit ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ ) ;
123	Double_t ProbGoodChiSqHit ( Double_t radius, Double_t searchRadiusRPhi, Double_t searchRadiusZ ) ;
124	Double_t ProbGoodChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ ) ;
125	Double_t ProbNullChiSqPlusConfHit ( Double_t radius, Double_t leff, Double_t searchRadiusRPhi, Double_t searchRadiusZ ) ;
126
127	// Howard W. hit distribution and convolution integral
128	Double_t Dist ( Double_t Z, Double_t radius ) ;
129	Double_t HitDensity ( Double_t radius ) ;
130	Double_t UpcHitDensity ( Double_t radius ) ;
131	Double_t IntegratedHitDensity ( Double_t multiplicity, Double_t radius ) ;
132	Double_t OneEventHitDensity ( Double_t multiplicity, Double_t radius ) const ;
133	Double_t D0IntegratedEfficiency( Double_t pt, Double_t corrEfficiency[][400] ) const ;
134
135	TGraph* GetGraphMomentumResolution(Int_t color, Int_t linewidth=1);
136	TGraph* GetGraphPointingResolution(Int_t axis,Int_t color, Int_t linewidth=1);
137	TGraph* GetGraphPointingResolutionTeleEqu(Int_t axis,Int_t color, Int_t linewidth=1);
138	TGraph* GetGraphLayerInfo(Int_t plot, Int_t color, Int_t linewidth=1);
139
140	TGraph* GetGraphImpactParam(Int_t mode, Int_t axis, Int_t color, Int_t linewidth=1);
141
142	TGraph* GetGraphRecoEfficiency(Int_t particle, Int_t color, Int_t linewidth=1);
143	TGraph* GetGraphRecoFakes(Int_t particle,Int_t color, Int_t linewidth);
144	TGraph* GetGraphRecoPurity(Int_t particle,Int_t color, Int_t linewidth);
145
146	TGraph* GetGraph(Int_t number, Int_t color, Int_t linewidth=1);
147
148	void MakeStandardPlots(Bool_t add =0, Int_t color=1, Int_t linewidth=1,Bool_t onlyPionEff=0);
149
150	// method to extend AliExternalTrackParam functionality
151	static Bool_t GetXatLabR(AliExternalTrackParam* tr,Double_t r,Double_t &x, Double_t bz, Int_t dir=0);
152	static Bool_t PropagateToR(AliExternalTrackParam* trc, double r, double b, int dir=0);
153	Double_t* PrepareEffFakeKombinations(TMatrixD probKomb, TMatrixD probLay);
154
155	protected:
156
157	Int_t fNumberOfLayers; // total number of layers in the model
158	Int_t fNumberOfActiveLayers; // number of active layers in the model
159	Int_t fNumberOfActiveITSLayers; // number of active ITS layers in the model
160	TList fLayers; // List of layer pointers
161	Float_t fBField; // Magnetic Field in Tesla
162	Float_t fLhcUPCscale; // UltraPeripheralElectrons: scale from RHIC to LHC
163	Float_t fIntegrationTime; // electronics integration time
164	Float_t fConfLevel; // Confidence Level for the tracking
165	Float_t fAvgRapidity; // rapidity of the track (= mean)
166	Float_t fParticleMass; // Particle used for tracking. Standard: mass of pion
167	Double_t fMaxRadiusSlowDet; // Maximum radius for slow detectors. Fast detectors
168	// and only fast detectors reside outside this radius.
169	Int_t fAtLeastCorr; // min. number of correct hits for the track to be "good"
170	Int_t fAtLeastFake; // min. number of fake hits for the track to be "fake"
171
172	Double_t fMaxSeedRadius; // max seeding radius, e.g. to exclude TRD
173	Double_t fptScale; // used for maxPt
174
175	Int_t fdNdEtaCent; // Multiplicity
176
177	enum {kMaxNumberOfDetectors = 200};
178
179	Double_t fTransMomenta[kNptBins]; // array of transverse momenta
180	Double_t fMomentumRes[kNptBins]; // array of momentum resolution
181	Double_t fResolutionRPhi[kNptBins]; // array of rphi resolution
182	Double_t fResolutionZ[kNptBins]; // array of z resolution
183	Double_t fDetPointRes[kMaxNumberOfDetectors][kNptBins]; // array of rphi resolution per layer
184	Double_t fDetPointZRes[kMaxNumberOfDetectors][kNptBins]; // array of z resolution per layer
185	Double_t fEfficiency[3][kNptBins]; // efficiency for different particles
186	Double_t fFake[3][kNptBins]; // fake prob for different particles
187
188	Int_t kDetLayer; // layer for which a few more details are extracted
189	Double_t fResolutionRPhiLay[kNptBins]; // array of rphi resolution
190	Double_t fResolutionZLay[kNptBins]; // array of z resolution
191	Double_t fEfficProlongLay[kNptBins]; // array of z resolution
192
193	ClassDef(DetectorK,1);
194	};
195
196	#endif
197