Online-Offline channel correspondence
[u/mrichter/AliRoot.git] / AD / ADbase / AliADDataDCS.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
17 #include "AliADDataDCS.h"
18
19 #include "AliDCSValue.h"
20 #include "AliLog.h"
21 #include "AliADConst.h"
22
23 #include <TGraph.h>
24 #include <TAxis.h>
25 #include <TCanvas.h>
26 #include <TTimeStamp.h>
27 #include <TMap.h>
28 #include <TString.h>
29 #include <TObjString.h>
30 #include <TH1F.h>
31
32 class TH2;
33 class AliCDBMetaData;
34 class TDatime;
35
36 // AliADDataDCS class
37 // main aim to introduce the aliases for the AD DCS
38 // data points to be then
39 // stored in the OCDB, and to process them. 
40 // ProcessData() method called by ADPreprocessor
41
42 ClassImp(AliADDataDCS)
43
44 //_____________________________________________________________________________
45 AliADDataDCS::AliADDataDCS():
46         TObject(),
47         fRun(0),
48         fStartTime(0),
49         fEndTime(0),
50         fDaqStartTime(0),
51         fDaqEndTime(0),
52         fCtpStartTime(0),
53         fCtpEndTime(0),
54     fGraphs("TGraph",kNGraphs),
55         fFEEParameters(NULL),
56         fIsProcessed(kFALSE)
57
58 {
59   // Default constructor
60         for(int i=0;i<kNHvChannel;i++) {
61                 fDeadChannel[i] = kFALSE;
62                 fMeanHV[i]      = 100.0;
63                 fWidthHV[i]     = 0.0; 
64                 fHv[i]          = NULL;
65         }
66 }
67
68 //_____________________________________________________________________________
69 AliADDataDCS::AliADDataDCS(Int_t nRun, UInt_t startTime, UInt_t endTime, UInt_t daqStartTime, UInt_t daqEndTime, UInt_t ctpStartTime, UInt_t ctpEndTime):
70         TObject(),
71         fRun(nRun),
72         fStartTime(startTime),
73         fEndTime(endTime),
74         fDaqStartTime(daqStartTime),
75         fDaqEndTime(daqEndTime),
76         fCtpStartTime(ctpStartTime),
77         fCtpEndTime(ctpEndTime),
78         fGraphs("TGraph",kNGraphs),
79         fFEEParameters(new TMap()),
80         fIsProcessed(kFALSE)
81
82 {
83
84   // constructor with arguments
85         for(int i=0;i<kNHvChannel;i++) {
86          fDeadChannel[i] = kFALSE;        
87          fMeanHV[i]      = 100.0;
88      fWidthHV[i]     = 0.0; 
89         }
90         AliInfo(Form("\n\tRun %d \n\tTime Created %s \n\tTime Completed %s \n\tDAQ start %s \n\tDAQ end %s \n\tCTP start %s \n\tCTP end %s   ", nRun,
91                 TTimeStamp(startTime).AsString(),
92                 TTimeStamp(endTime).AsString(),
93                 TTimeStamp(daqStartTime).AsString(),
94                 TTimeStamp(daqEndTime).AsString(),
95                 TTimeStamp(ctpStartTime).AsString(),
96                 TTimeStamp(ctpEndTime).AsString()
97                 ));
98         
99         fFEEParameters->SetOwnerValue();
100         Init();
101
102 }
103
104 //_____________________________________________________________________________
105 AliADDataDCS::~AliADDataDCS() {
106
107   // destructor
108   fGraphs.Clear("C");
109   delete fFEEParameters;
110
111 }
112
113 //_____________________________________________________________________________
114 Bool_t AliADDataDCS::ProcessData(TMap& aliasMap){
115
116   // method to process the data
117   Bool_t success = kTRUE;
118
119   if(!(fAliasNames[0])) Init();
120
121   TObjArray *aliasArr;
122   AliDCSValue* aValue;
123
124   // starting loop on aliases
125   for(int iAlias=0; iAlias<kNAliases; iAlias++){
126
127     aliasArr = (TObjArray*) aliasMap.GetValue(fAliasNames[iAlias].Data());
128     if(!aliasArr){
129       AliError(Form("Alias %s not found!", fAliasNames[iAlias].Data()));
130       success = kFALSE;
131       continue;
132     }
133
134     //Introduce(iAlias, aliasArr);
135     
136     if(aliasArr->GetEntries()<2){
137       AliWarning(Form("Alias %s has just %d entries!",
138                     fAliasNames[iAlias].Data(),aliasArr->GetEntries()));
139     }
140     
141     TIter iterarray(aliasArr);
142         
143     if(iAlias<kNHvChannel){ // Treating HV values
144         Int_t nentries = aliasArr->GetEntries();
145
146         Double_t *times = new Double_t[nentries];
147         Double_t *values = new Double_t[nentries];
148
149         UInt_t iValue=0;
150         Float_t variation = 0.0;
151
152         while((aValue = (AliDCSValue*) iterarray.Next())) {
153                         UInt_t currentTime = aValue->GetTimeStamp();
154                         //if(currentTime>fCtpEndTime) break; //What is this for?
155
156                         values[iValue] = aValue->GetFloat();
157                         times[iValue] = (Double_t) (currentTime);
158                         
159                         if(iValue>0) {
160                                 if(values[iValue-1]>0.) variation = TMath::Abs(values[iValue]-values[iValue-1])/values[iValue-1];
161                                 if(variation > 0.01) fDeadChannel[kOfflineChannel[iAlias]] = kTRUE;
162                         }
163                         fHv[iAlias]->Fill(values[iValue]);
164                         printf("%s : %s : %f Dead=%d\n",fAliasNames[iAlias].Data(),TTimeStamp(currentTime).AsString(),values[iValue],fDeadChannel[kOfflineChannel[iAlias]]);
165                         iValue++;
166         }      
167         CreateGraph(iAlias, aliasArr->GetEntries(), times, values); // fill graphs 
168
169         // calculate mean and rms of the first two histos
170         // and convert index to aliroot channel
171         Int_t iChannel     = kOfflineChannel[iAlias];   
172         fMeanHV[iChannel]  = fHv[iAlias]->GetMean();
173         fWidthHV[iChannel] = fHv[iAlias]->GetRMS();
174
175         delete[] values;
176         delete[] times; 
177         } else { // Treating FEE Parameters
178                 AliDCSValue * lastVal = NULL;
179                 while((aValue = (AliDCSValue*) iterarray.Next())) lastVal = aValue; // Take only the last value
180                 fFEEParameters->Add(new TObjString(fAliasNames[iAlias].Data()),lastVal);
181         }      
182   }
183   
184   fIsProcessed=kTRUE;
185
186   return success;
187 }
188
189 //_____________________________________________________________________________
190 void AliADDataDCS::Init(){
191
192   // initialization of aliases and DCS data
193
194   TString sindex;
195   int iAlias = 0;
196   
197   for(int iPM = 0; iPM<16 ; iPM++){
198                 fAliasNames[iAlias] = Form("AD0/HV/PM%d",iPM);
199                         
200                 fHv[iAlias] = new TH1F(fAliasNames[iAlias].Data(),fAliasNames[iAlias].Data(), 3000, kHvMin, kHvMax);
201                 fHv[iAlias]->GetXaxis()->SetTitle("Hv");
202                 iAlias++;
203   }
204   // CCIU Parameters
205         
206   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BBAThreshold";
207   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BBCThreshold";
208   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BGAThreshold";
209   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BGCThreshold";
210   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BBAForBGThreshold";
211   fAliasNames[iAlias++] = "AD0/FEE/CCIU/BBCForBGThreshold";
212   fAliasNames[iAlias++] = "AD0/FEE/CCIU/MultADAThrLow";
213   fAliasNames[iAlias++] = "AD0/FEE/CCIU/MultADAThrHigh";
214   fAliasNames[iAlias++] = "AD0/FEE/CCIU/MultADCThrLow";
215   fAliasNames[iAlias++] = "AD0/FEE/CCIU/MultADCThrHigh";
216   for(int i=1;i<=5;i++) {
217         fAliasNames[iAlias] = "AD0/FEE/CCIU/TriggerSelect";
218         fAliasNames[iAlias++] += Form("%d",i);
219   }
220   //Trigger parameters
221   for(int iCIU = 0; iCIU<2 ; iCIU++){
222                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/Clk1Win1",iCIU);
223                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/Clk1Win2",iCIU);
224                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/Clk2Win1",iCIU);
225                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/Clk2Win2",iCIU);
226                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DelayClk1Win1",iCIU);
227                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DelayClk1Win2",iCIU);
228                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DelayClk2Win1",iCIU);
229                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DelayClk2Win2",iCIU);
230                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/LatchWin1",iCIU);
231                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/LatchWin2",iCIU);
232                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/ResetWin1",iCIU);
233                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/ResetWin2",iCIU);
234                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/PedestalSubtraction",iCIU);
235   }
236   for(int iCIU = 0; iCIU<2 ; iCIU++){
237           for(int iCh=1;iCh<=8;iCh++){
238                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/EnableCharge%d",iCIU,iCh);
239                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/EnableTiming%d",iCIU,iCh);
240                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/PedEven%d",iCIU,iCh);
241                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/PedOdd%d",iCIU,iCh);
242                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/PedCutEven%d",iCIU,iCh);
243                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/PedCutOdd%d",iCIU,iCh);
244           }
245   }
246  // Time Resolution Parameters  
247   for(int iCIU = 0; iCIU<2 ; iCIU++){
248                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/TimeResolution",iCIU);
249                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/WidthResolution",iCIU);
250   }
251   // HPTDC parameters
252   for(int iCIU = 0; iCIU<2 ; iCIU++){
253                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/MatchWindow",iCIU);
254                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/SearchWindow",iCIU);
255                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/TriggerCountOffset",iCIU);
256                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/RollOver",iCIU);
257   }
258
259   for(int iCIU = 0; iCIU<2 ; iCIU++){
260           for(int iCh=1;iCh<=8;iCh++){
261                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DelayHit%d",iCIU,iCh);
262           }
263   }
264
265   for(int iCIU = 0; iCIU<2 ; iCIU++){
266           for(int iCh=1;iCh<=8;iCh++){
267                 fAliasNames[iAlias++] = Form("AD0/FEE/CIU%d/DiscriThr%d",iCIU,iCh);
268           }
269   }
270
271   if(iAlias!=kNAliases) 
272               AliError(Form("Number of DCS Aliases defined not correct"));
273
274 }
275
276 //_____________________________________________________________________________
277 void AliADDataDCS::Introduce(UInt_t numAlias, const TObjArray* aliasArr)const
278 {
279
280   // method to introduce new aliases
281
282   int entries=aliasArr->GetEntries();
283   AliInfo(Form("************ Alias: %s **********",fAliasNames[numAlias].Data()));
284   AliInfo(Form("        %d DP values collected",entries));
285
286 }
287
288 //_____________________________________________________________________________
289 void AliADDataDCS::CreateGraph(int i, int dim, const Double_t *x, const Double_t *y)
290 {
291
292    // Create graphics
293    
294    TGraph *gr = new(fGraphs[fGraphs.GetEntriesFast()]) TGraph(dim, x, y);
295
296    gr->GetXaxis()->SetTimeDisplay(1);
297    gr->SetTitle(fAliasNames[i].Data());
298
299    AliInfo(Form("Array entries: %d",fGraphs.GetEntriesFast()));
300
301 }
302
303
304 //_____________________________________________________________________________
305 void AliADDataDCS::Draw(const Option_t* /*option*/)
306 {
307 // Draw all histos and graphs
308
309   if(!fIsProcessed) return;
310
311   if(fGraphs.GetEntries()==0)  return;
312   
313   TCanvas *cHV = new TCanvas("AD0_HV","AD0_HV");
314   cHV->Divide(4,4);
315   
316   for(int iPM = 0; iPM<16 ; iPM++){
317         cHV->cd(iPM+1);
318         ((TGraph*) fGraphs.UncheckedAt(iPM))->SetMarkerStyle(20);
319         ((TGraph*) fGraphs.UncheckedAt(iPM))->Draw("ALP");                      
320         }
321 }
322