]>
Commit | Line | Data |
---|---|---|
0ccdab7b | 1 | #ifndef ALIFMDESDFIXER_H |
2 | #define ALIFMDESDFIXER_H | |
3 | #include <TObject.h> | |
4 | #include <TBits.h> | |
5 | class TH1; | |
6 | class TList; | |
7 | class AliESDFMD; | |
8 | ||
9 | /** | |
10 | * Class to fix up an ESD object for various small issues. | |
11 | * | |
12 | * @par Input: | |
13 | * - AliESDFMD object - from reconstruction | |
14 | * - Ip z coordinate - from reconstruction | |
15 | * - Reco noise factor - Assumed noise factor used in reconstruction | |
16 | * | |
17 | * @par Output: | |
18 | * - The same AliESDFMD object but modified | |
19 | * | |
20 | * @par Corrections used | |
21 | * - AliFMDCorrNoiseGain if the noise correction is enabled | |
22 | * | |
23 | * @par Histograms | |
24 | * - Change in @f$\Delta@f$ due to noise correction | |
25 | * - Number of channels declared dead | |
26 | * - Change in @f$\eta@f$ | |
27 | * | |
28 | * @ingroup pwglf_forward_algo | |
29 | * @ingroup pwglf_forward_aod | |
30 | */ | |
31 | class AliFMDESDFixer : public TObject | |
32 | { | |
33 | public: | |
34 | /** | |
35 | * Default CTOR - for ROOT I/O - do not use | |
36 | */ | |
37 | AliFMDESDFixer(); | |
38 | ||
39 | /** | |
40 | * User CTOR | |
41 | * | |
42 | * @param name Dummy argument | |
43 | */ | |
44 | AliFMDESDFixer(const char* name); | |
45 | /** | |
46 | * Get name of object | |
47 | * | |
48 | * @return Always the same | |
49 | */ | |
50 | const char* GetName() const { return "fmdESDFixer"; } | |
51 | /** | |
52 | * Create output objects | |
53 | */ | |
54 | void CreateOutputObjects(TList* l); | |
55 | /** | |
56 | * Fix the ESD object | |
57 | * | |
58 | * @param esd ESD object | |
59 | * @param ipZ IP Z--coordinate | |
60 | */ | |
61 | void Fix(AliESDFMD& esd, Double_t ipZ); | |
62 | ||
63 | /** | |
64 | * @{ | |
65 | * @name Noise suppression fix | |
66 | */ | |
67 | /** | |
68 | * Set the factor assumed to be used in the reconstruction. If this | |
69 | * is set way high (>=4) then this corrector will effectively be | |
70 | * disabled. | |
71 | * | |
72 | * @param f | |
73 | */ | |
74 | void SetRecoNoiseFactor(Int_t f) { fRecoFactor = f; } | |
75 | void SetMaxNoiseCorrection(Double_t x) { fMaxNoiseCorr = x; } | |
76 | /** | |
77 | * Get the factor assumed in reconstruction pass | |
78 | * | |
79 | * @return factor assumed in reconstruction pass | |
80 | */ | |
81 | Int_t GetRecoNoiseFactor() const { return fRecoFactor; } | |
82 | /** | |
83 | * Check if we're using the noise correction. | |
84 | * | |
85 | * @return true if fRecoFactor < 4 | |
86 | */ | |
87 | Bool_t IsUseNoiseCorrection() const { return fRecoFactor < 4; } | |
88 | /** | |
89 | * Find the target noise factor | |
90 | * | |
91 | * @param esd ESD object. | |
92 | * @param check If true, also check for correction object | |
93 | * | |
94 | * @return Needed noise factor, or 0 or less if no correction is needed | |
95 | */ | |
96 | Int_t FindTargetNoiseFactor(const AliESDFMD& esd, Bool_t check=true) const; | |
97 | /* @} */ | |
98 | ||
99 | /** | |
100 | * @{ | |
101 | * @name Recalculations | |
102 | */ | |
103 | /** | |
104 | * In case of a displaced vertices recalculate @f$\eta@f$ and angle | |
105 | * correction | |
106 | * | |
107 | * @param use recalculate or not | |
108 | * | |
109 | */ | |
110 | void SetRecalculateEta(Bool_t use) { fRecalculateEta = use; } | |
111 | /* @} */ | |
112 | ||
113 | /** | |
114 | * @{ | |
115 | * @name Special treatmeant of invalid signals | |
116 | */ | |
117 | /** | |
118 | * Set whether to consider invalid multiplicities as null (or empty) | |
119 | * signal. | |
120 | * | |
121 | * @param flag If true, count invalids as empty | |
122 | */ | |
123 | void SetInvalidIsEmpty(Bool_t flag) { fInvalidIsEmpty = flag; } | |
124 | /** @} */ | |
125 | ||
126 | /** | |
127 | * @{ | |
128 | * @name Dead strip handling | |
129 | */ | |
130 | /** | |
131 | * Add a dead strip | |
132 | * | |
133 | * @param d Detector | |
134 | * @param r Ring | |
135 | * @param s Sector | |
136 | * @param t Strip | |
137 | */ | |
138 | void AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t); | |
139 | /** | |
140 | * Add a dead region in a detector ring | |
141 | * | |
142 | * @param d Detector | |
143 | * @param r Ring | |
144 | * @param s1 First sector (inclusive) | |
145 | * @param s2 Last sector (inclusive) | |
146 | * @param t1 First strip (inclusive) | |
147 | * @param t2 Last strip (inclusive) | |
148 | */ | |
149 | void AddDeadRegion(UShort_t d, Char_t r, UShort_t s1, UShort_t s2, | |
150 | UShort_t t1, UShort_t t2); | |
151 | /** | |
152 | * Add dead strips from a script. The script is supposed to accept | |
153 | * a pointer to this object (AliFMDSharingFilter) and then call | |
154 | * AddDead or AddDeadRegion as needed. | |
155 | * | |
156 | * @code | |
157 | * void deadstrips(AliFMDSharingFilter* filter) | |
158 | * { | |
159 | * filter->AddDead(...); | |
160 | * // ... and so on | |
161 | * } | |
162 | * @endcode | |
163 | * | |
164 | * @param script The script to read dead strips from. | |
165 | */ | |
166 | void AddDead(const Char_t* script); | |
167 | /* @} */ | |
168 | ||
169 | void Print(Option_t* option="") const; | |
170 | protected: | |
171 | AliFMDESDFixer(const AliFMDESDFixer&); | |
172 | AliFMDESDFixer& operator=(const AliFMDESDFixer&); | |
173 | ||
174 | /** | |
175 | * @{ | |
176 | * @name Worker functions | |
177 | */ | |
178 | /** | |
179 | * Check if a strip is marked as dead | |
180 | * | |
181 | * @param d Detector | |
182 | * @param r Ring | |
183 | * @param s Sector | |
184 | * @param t Strip | |
185 | * | |
186 | * @return true if dead | |
187 | */ | |
188 | virtual Bool_t IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const; | |
189 | /** | |
190 | * Possibly raise a strip from the dead or kill it | |
191 | * | |
192 | * @param d Detector | |
193 | * @param r Ring | |
194 | * @param s Sector | |
195 | * @param t Strip | |
196 | * @param m Multiplicity | |
197 | * | |
198 | * @return true if this was killed | |
199 | */ | |
200 | Bool_t CheckDead(UShort_t d, Char_t r, UShort_t s, UShort_t t, Double_t& m); | |
201 | /** | |
202 | * Re-calculate eta and correct the multiplicity accordingly | |
203 | * | |
204 | * @param d Detector | |
205 | * @param r Ring | |
206 | * @param s Sector | |
207 | * @param t Strip | |
208 | * @param vz Ip Z-coordinate | |
209 | * @param mult In/out multiplicity | |
210 | * @param eta In/out eta | |
211 | * @param cosTheta On return, the cosine of theta or null | |
212 | */ | |
213 | void RecalculateEta(UShort_t d, Char_t r, UShort_t s, UShort_t t, | |
214 | Double_t vz, Double_t& mult, Double_t& eta, | |
215 | Double_t& cosTheta); | |
216 | /** | |
217 | * Correct for noise suppression | |
218 | * | |
219 | * @param f Factor to apply | |
220 | * @param c Correction | |
221 | * @param cosTheta Cosine to theta | |
222 | * @param mult In/Out multiplity | |
223 | * | |
224 | * @return true if signal is good, otherwise false | |
225 | */ | |
226 | Bool_t NoiseCorrect(Int_t f, Double_t c, Double_t cosTheta, Double_t& mult); | |
227 | ||
228 | Int_t fRecoFactor; // Noise factor used in Reco | |
229 | Double_t fMaxNoiseCorr; // If noise corr above this, flag as dead | |
230 | Bool_t fRecalculateEta; // Whether to recalc eta and angle cor (disp vtx) | |
231 | TBits fXtraDead; // List of extra dead channels | |
232 | Bool_t fHasXtraDead; // Whether we have xtra dead channels | |
233 | Bool_t fInvalidIsEmpty; // Consider kInvalidMult as zero | |
234 | TH1* fNoiseChange; // Diagnostics | |
235 | TH1* fEtaChange; // Diagnostics | |
236 | TH1* fDeadChange; // Diagnostics | |
237 | ||
238 | ClassDef(AliFMDESDFixer,1); // Fix FMD ESD object for issues | |
239 | }; | |
240 | ||
241 | #endif | |
242 | // | |
243 | // EOF | |
244 | // | |
245 |