Using detector quality flag (taken from ALICE logbook) to decide whether to rpodcue...
[u/mrichter/AliRoot.git] / VZERO / AliVZEROTriggerMask.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17 //-------------------------------------------
18 // Class : AliVZEROTriggerMask
19 //
20 // Fill up the trigger mask word.
21 //
22
23 #include <Riostream.h>
24 #include <TGeoManager.h>
25 #include <TGeoMatrix.h>
26 #include <TGeoPhysicalNode.h>
27 #include <TF1.h>
28 #include <TMath.h>
29
30 #include <AliGeomManager.h>
31 #include "AliLog.h"
32 #include "AliVZEROTriggerMask.h"
33 #include "AliVZEROConst.h"
34 #include "AliVZEROCalibData.h"
35 #include "AliESDVZERO.h"
36 #include "AliVZEROReconstructor.h"
37
38 //______________________________________________________________________
39 ClassImp(AliVZEROTriggerMask)
40
41 //______________________________________________________________________
42
43 AliVZEROTriggerMask::AliVZEROTriggerMask()
44   :TObject(),
45    fV0ADist(0),
46    fV0CDist(0),
47    fRecoParam(NULL)
48 {
49   // Default constructor
50   //
51   Float_t zV0A = TMath::Abs(GetZPosition("VZERO/V0A"));
52   Float_t zV0C = TMath::Abs(GetZPosition("VZERO/V0C"));
53
54   // distance in time units from nominal vertex to V0
55   fV0ADist = zV0A/TMath::Ccgs()*1e9;
56   fV0CDist = zV0C/TMath::Ccgs()*1e9;
57 }
58
59 //________________________________________________________________________________
60 Double_t AliVZEROTriggerMask::GetZPosition(const char* symname){
61 // Get the global z coordinate of the given V0 alignable volume
62 //
63   Double_t *tr;
64   TGeoPNEntry *pne = gGeoManager->GetAlignableEntry(symname);
65   if (!pne) {
66     AliFatalClass(Form("TGeoPNEntry with symbolic name %s does not exist!",symname));
67     return 0;
68   }
69
70   TGeoPhysicalNode *pnode = pne->GetPhysicalNode();
71   if(pnode){
72           TGeoHMatrix* hm = pnode->GetMatrix();
73            tr = hm->GetTranslation();
74   }else{
75           const char* path = pne->GetTitle();
76           if(!gGeoManager->cd(path)){
77                   AliFatalClass(Form("Volume path %s not valid!",path));
78                   return 0;
79           }
80          tr = gGeoManager->GetCurrentMatrix()->GetTranslation();
81   }
82   return tr[2];
83
84 }
85
86 //________________________________________________________________________________
87
88
89 void AliVZEROTriggerMask::FillMasks(AliESDVZERO *esdV0,
90                                     AliVZEROCalibData *cal,
91                                     TF1 *slewing)
92 {
93   // Fill up the trigger mask word
94   // using the TDC data (already corrected for
95   // slewing and misalignment between channels)
96
97   esdV0->SetBit(AliESDVZERO::kTriggerBitsFilled,kTRUE);
98   esdV0->SetBit(AliESDVZERO::kDecisionFilled,kTRUE);
99
100   const Float_t p1 = 2.50; // photostatistics term in the time resolution
101   const Float_t p2 = 3.00; // slewing related term in the time resolution
102
103   // loop over vzero channels
104   Float_t timeAW = 0,timeCW = 0;
105   Float_t weightA = 0,weightC = 0;
106   Int_t ntimeA = 0, ntimeC = 0;
107   Double_t timesA[32], timesC[32];
108   Double_t wA[32],wC[32];
109   Int_t indA[32], indC[32];
110   for (Int_t i = 0; i < 64; ++i) {
111     Float_t adc = esdV0->GetAdc(i);
112     if (adc > GetRecoParam()->GetAdcThresHold()) {
113       Float_t tdc = esdV0->GetTime(i);
114       if (tdc > (AliVZEROReconstructor::kInvalidTime + 1e-6)) {
115         Float_t nphe = adc*kChargePerADC/(cal->GetGain(i)*TMath::Qe());
116         Float_t timeErr = TMath::Sqrt(kIntTimeRes*kIntTimeRes+
117                                       p1*p1/nphe+
118                                       p2*p2*(slewing->GetParameter(0)*slewing->GetParameter(1))*(slewing->GetParameter(0)*slewing->GetParameter(1))*
119                                       TMath::Power(adc/cal->GetCalibDiscriThr(i,kTRUE),2.*(slewing->GetParameter(1)-1.))/
120                                       (cal->GetCalibDiscriThr(i,kTRUE)*cal->GetCalibDiscriThr(i,kTRUE)));
121
122         if (i < 32) { // in V0C
123           timesC[ntimeC] = tdc;
124           wC[ntimeC] = 1.0/(timeErr*timeErr);
125           indC[ntimeC] = i;
126           ntimeC++;
127           timeCW += tdc/(timeErr*timeErr);
128           weightC += 1.0/(timeErr*timeErr);
129         }
130         else { // in V0A
131           timesA[ntimeA] = tdc;
132           wA[ntimeA] = 1.0/(timeErr*timeErr);
133           indA[ntimeA] = i - 32;
134           ntimeA++;
135           timeAW += tdc/(timeErr*timeErr);
136           weightA += 1.0/(timeErr*timeErr);
137         }
138       }
139     }
140   } // end of loop over channels
141
142   if (weightA > 0) timeAW = timeAW/weightA;
143   else timeAW = AliVZEROReconstructor::kInvalidTime;
144
145   if (weightC > 0) timeCW = timeCW/weightC;
146   else timeCW = AliVZEROReconstructor::kInvalidTime;
147
148   esdV0->SetBit(AliESDVZERO::kRobustMeanTime,kTRUE);
149
150   Double_t medianTimeA = AliVZEROReconstructor::kInvalidTime;
151   if (ntimeA > 0) medianTimeA = TMath::Median(ntimeA,timesA,wA);
152   Double_t medianTimeC = AliVZEROReconstructor::kInvalidTime;
153   if (ntimeC > 0) medianTimeC = TMath::Median(ntimeC,timesC,wC);
154
155   Float_t robTimeAW = 0,robTimeCW = 0;
156   Float_t robWeightA = 0,robWeightC = 0;
157   Int_t nrobTimeA = 0, nrobTimeC = 0;
158   Int_t robIndA[32], robIndC[32];
159   for(Int_t i = 0; i < ntimeA; ++i) {
160     AliDebug(1,Form("ChannelsAResiduals %f %f %d",timesA[i]-medianTimeA,1./TMath::Sqrt(wA[i]),ntimeA));
161     if (TMath::Abs(timesA[i]-medianTimeA) < GetRecoParam()->GetMaxResid()/TMath::Sqrt(wA[i])) {
162       robIndA[nrobTimeA] = indA[i];
163       nrobTimeA++;
164       robTimeAW += timesA[i]*wA[i];
165       robWeightA += wA[i];
166     }
167   }
168   for(Int_t i = 0; i < ntimeC; ++i) {
169     AliDebug(1,Form("ChannelsCResiduals %f %f %d",timesC[i]-medianTimeC,1./TMath::Sqrt(wC[i]),ntimeC));
170     if (TMath::Abs(timesC[i]-medianTimeC) < GetRecoParam()->GetMaxResid()/TMath::Sqrt(wC[i])) {
171       robIndC[nrobTimeC] = indC[i];
172       nrobTimeC++;
173       robTimeCW += timesC[i]*wC[i];
174       robWeightC += wC[i];
175     }
176   }
177
178   if (robWeightA > 0) robTimeAW = robTimeAW/robWeightA;
179   else robTimeAW = AliVZEROReconstructor::kInvalidTime;
180
181   if (robWeightC > 0) robTimeCW = robTimeCW/robWeightC;
182   else robTimeCW = AliVZEROReconstructor::kInvalidTime;
183
184   AliDebug(1,Form("V0timesA %f %f %f %f %d",timeAW,(weightA > 0) ? (1./TMath::Sqrt(weightA)) : 0,
185                   medianTimeA,robTimeAW,ntimeA));
186   AliDebug(1,Form("V0timesC %f %f %f %f %d",timeCW,(weightC > 0) ? (1./TMath::Sqrt(weightC)) : 0,
187                   medianTimeC,robTimeCW,ntimeC));
188
189   esdV0->SetV0ATime(robTimeAW);
190   esdV0->SetV0CTime(robTimeCW);
191   esdV0->SetV0ATimeError((robWeightA > 0) ? (1./TMath::Sqrt(robWeightA)) : 0);
192   esdV0->SetV0CTimeError((robWeightC > 0) ? (1./TMath::Sqrt(robWeightC)) : 0);
193
194   esdV0->SetV0ADecision(AliESDVZERO::kV0Empty);
195   esdV0->SetV0CDecision(AliESDVZERO::kV0Empty);
196
197   if (robTimeAW > (fV0ADist + GetRecoParam()->GetTimeWindowBBALow()) &&
198       robTimeAW < (fV0ADist + GetRecoParam()->GetTimeWindowBBAUp())) 
199     esdV0->SetV0ADecision(AliESDVZERO::kV0BB);
200   else if (robTimeAW > (-fV0ADist + GetRecoParam()->GetTimeWindowBGALow()) &&
201            robTimeAW < (-fV0ADist + GetRecoParam()->GetTimeWindowBGAUp()))
202     esdV0->SetV0ADecision(AliESDVZERO::kV0BG);
203   else if (robTimeAW > (AliVZEROReconstructor::kInvalidTime + 1e-6))
204     esdV0->SetV0ADecision(AliESDVZERO::kV0Fake);
205
206   if (robTimeCW > (fV0CDist + GetRecoParam()->GetTimeWindowBBCLow()) &&
207       robTimeCW < (fV0CDist + GetRecoParam()->GetTimeWindowBBCUp())) 
208     esdV0->SetV0CDecision(AliESDVZERO::kV0BB);
209   else if (robTimeCW > (-fV0CDist + GetRecoParam()->GetTimeWindowBGCLow()) &&
210            robTimeCW < (-fV0CDist + GetRecoParam()->GetTimeWindowBGCUp()))
211     esdV0->SetV0CDecision(AliESDVZERO::kV0BG);
212   else if (robTimeCW > (AliVZEROReconstructor::kInvalidTime + 1e-6))
213     esdV0->SetV0CDecision(AliESDVZERO::kV0Fake);
214
215   UInt_t aBBtriggerV0A = 0; // bit mask for Beam-Beam trigger in V0A
216   UInt_t aBGtriggerV0A = 0; // bit mask for Beam-Gas trigger in V0A
217   UInt_t aBBtriggerV0C = 0; // bit mask for Beam-Beam trigger in V0C
218   UInt_t aBGtriggerV0C = 0; // bit mask for Beam-Gas trigger in V0C
219
220   for(Int_t i = 0; i < nrobTimeA; ++i) {
221     if (esdV0->GetV0ADecision() == AliESDVZERO::kV0BB)
222       aBBtriggerV0A |= (1 << (robIndA[i]));
223     else if (esdV0->GetV0ADecision() == AliESDVZERO::kV0BG)
224       aBGtriggerV0A |= (1 << (robIndA[i]));
225   }
226
227   for(Int_t i = 0; i < nrobTimeC; ++i) {
228     if (esdV0->GetV0CDecision() == AliESDVZERO::kV0BB)
229       aBBtriggerV0C |= (1 << (robIndC[i]));
230     else if (esdV0->GetV0CDecision() == AliESDVZERO::kV0BG)
231       aBGtriggerV0C |= (1 << (robIndC[i]));
232   }
233
234   esdV0->SetBBtriggerV0A(aBBtriggerV0A);
235   esdV0->SetBGtriggerV0A(aBGtriggerV0A);
236   esdV0->SetBBtriggerV0C(aBBtriggerV0C);
237   esdV0->SetBGtriggerV0C(aBGtriggerV0C);
238 }